Observational constraints on bimetric gravity
PPrepared for submission to JCAP
Observational constraints on bimetricgravity
Marcus H¨og˚as, Edvard M¨ortsell
The Oskar Klein Centre, Department of Physics, Stockholm University, SE 106 91, Stock-holm, SwedenE-mail: [email protected], [email protected]
Abstract.
Ghost-free bimetric gravity is a theory of two interacting spin-2 fields, one masslessand one massive, in addition to the standard matter particles and fields, thereby generalizingEinstein’s theory of general relativity. To parameterize the theory, we use five observableswith specific physical interpretations. We present, for the first time, observational constraintson these parameters that: (i) apply to the full theory, (ii) are consistent with a workingscreening mechanism (i.e., restoring general relativity locally), (iii) exhibit a continuous, real-valued background cosmology (without the Higuchi ghost). For the cosmological constraints,we use data sets from the cosmic microwave background, baryon acoustic oscillations, andtype Ia supernovae. Bimetric cosmology provides a good fit to data even for large valuesof the mixing angle between the massless and massive gravitons. Interestingly, the best-fitmodel is a self-accelerating solution where the accelerated expansion is due to the dynamicalmassive spin-2 field, without a cosmological constant. Due to the screening mechanism, themodels are consistent with local tests of gravity such as solar system tests and gravitationallensing by galaxies. We also comment on the possibility of alleviating the Hubble tensionwith this theory. a r X i v : . [ g r- q c ] J a n ontents D.1 Numerical errors 27D.2 Scanning details 29
Ghost-free bimetric gravity is a theory of two interacting spin-2 field, one massive and onemassless, in addition to the standard matter particles and fields. To achieve this setup, onemust introduce an auxiliary symmetric spin-2 field (i.e., a metric) f µν besides the physicalmetric g µν . The physical metric is the one to which standard matter particles and fieldscouple.Among the virtues of bimetric gravity are self-accelerating cosmological models, thatis without a cosmological constant [1–8]. For these models, the accelerated expansion ofthe Universe is due to the interaction between the spin-2 fields which contributes an extraterm, Ω DE , in the cosmological equations of motion. Interestingly, for these models, a smallvalue of Ω DE is technically natural in the sense of ’t Hooft, which means that its value isprotected from quantum corrections [9]. Bimetric cosmology can push the Hubble constantinferred from cosmic microwave background observations in the right direction in order to– 1 –ase the tension with local measurements [10]. The theory has a screening mechanism thatcan restore general relativity on solar system scales [11–13]. Among the challenges is theexistence of a gradient instability for linear perturbations around cosmological backgroundsolutions [14–29] and finding a stable (well-posed) form of the equations of motion in orderto obtain long-term numerical evolution of generic systems [30–37].Usually, bimetric theory is parameterized with five parameters, β n with n = 0 , , ..., β -parameters are not observables and the reportedconstraints on them depend on which convention is used, being different in different papers.Following Ref. [38], which generalizes the works of Refs. [29, 39], we use a parameterizationin terms of five observables: the mixing angle θ between the massless and massive gravi-tons, the graviton (Fierz–Pauli) mass m FP , the effective cosmological constant Ω Λ , and twoparameters, α and β (not to be confused with the β -parameters β n ), that appear in thescreening mechanism which cancels the extra gravitational forces on local scales (e.g., in thesolar system).Earlier results have shown the theory to be compatible with cosmological data from thecosmic microwave background (CMB), baryon acoustic oscillations (BAO), and supernovaetype Ia (SNIa) [2, 6, 8, 39–41], gravitational wave observations [42–45], solar system tests[13, 46], velocity dispersion and strong gravitational lensing in galaxies [12, 13, 47, 48], aswell as gravitational lensing by galaxy clusters [13]. Typically, in these works, only a subsetof observations has been used or the results apply only to a restricted set of models, that is,setting one or several of the β -parameters to zero. Also, the chosen set of parameters aretypically not compatible with a working screening mechanism or a continuous, real-valuedcosmology (without the Higuchi ghost). As shown in Ref. [38], these requirements can beformulated as a set of analytical constraints on the physical parameters. In this paper, forthe first time, we derive the observational constraints on the full theory while at the sametime guaranteeing that the above analytical constraints are satisfied. Notation.
For the most part, we use geometrized units where Newton’s gravitational con-stant and the speed of light are set to one, G = c = 1. Quantities constructed from theauxiliary metric f µν are denoted with tildes, otherwise constructed from the physical metric g µν . We use redshift z as a time variable, defined via the equation 1 + z = a /a where a isthe scale factor today. The equations of motion are constructed to avoid the Boulware–Deser ghost which plaguesgeneral theories of massive gravity [49, 50]. Assuming that there is only one matter sector,coupled to g µν , the equations read, G µν = κ g ( T µν + V µν ) , (cid:101) G µν = κ f (cid:101) V µν , (2.1)where G µν and (cid:101) G µν are the Einstein tensors, T µν is the standard matter stress–energy, and κ g and κ f are the gravitational constants of g µν and f µν , respectively. We denote the ratioof the gravitational constants, κ ≡ κ g /κ f . (2.2)The bimetric stress–energies V µν and (cid:101) V µν contain the metrics g µν and f µν as well as fiveconstant parameters β , ...β with dimension of curvature (see e.g. [38] for the mathematical– 2 –orm of V µν and (cid:101) V µν ). The bimetric stress–energies are conserved, ∇ µ V µν = (cid:101) ∇ µ (cid:101) V µν =0, which follows from conservation of matter stress–energy, ∇ µ T µν = 0, and the Bianchiidentities, ∇ µ G µν = (cid:101) ∇ µ (cid:101) G µν = 0. The β -parameters β , ..., β can be rescaled without affecting the physics. Hence, they arenot observables and the reported constraints on them depend on the choice of scaling, beingdifferent in different papers. To circumvent this problem, a subset of the physical parameterswas first introduced in Ref. [39] and then generalized to the full theory in Ref. [38]. Theseparameters are independent of the rescaling and are thus observable quantities. The frame-work applies to space-times where the metrics are proportional asymptotically (in space ortime), which we assume.At radial infinity of static, spherically symmetric solutions (i.e., local solutions, applica-ble on for example solar system scales) the metrics are proportional, that is g µν = c f µν with c = const. This is also true in the infinite future for the cosmological background solutions.We introduce dimensionless, rescaling invariant parameters [38], B n ≡ κ g β n c n /H , (3.1)expressing the β -parameters in units of the curvature scale defined by the Hubble constant.This is especially convenient for cosmological application since we expect β n ∼ H , hence B n ∼
1, if the theory is to exhibit novel features on cosmological scales. We can express thedimensionless physical parameters in terms of the B -parameters,tan θ = B + 3 B + 3 B + B B + 3 B + 3 B + B , (3.2a) m = ( B + 2 B + B ) / sin θ, (3.2b)Ω Λ = B B + B + B , (3.2c) α = − B + B B + 2 B + B , (3.2d) β = B B + 2 B + B . (3.2e)Inverting the relations, B = 3Ω Λ − sin θ m (3 + 3 α + β ) , (3.3a) B = sin θ m (1 + 2 α + β ) , (3.3b) B = − sin θ m ( α + β ) , (3.3c) B = sin θ m β, (3.3d) B = 3 tan θ Ω Λ + sin θ m ( − α − β ) . (3.3e)The physical parameters have specific physical interpretations: θ ∈ [0 , π/
2] is the mixingangle between the mass eigenstates and the metrics (cf. the mixing between neutrino massand flavor eigenstates), m FP > H = 100 h km / s / Mpc = 2 . h × − eV /c , Ω Λ > igure 1 : Exclusion plot in the αβ -plane from demanding a working screening mechanism(green), a continuous background cosmology without the Higuchi ghost (blue), and a positivematter density (red). Here, θ (cid:39) ◦ . The dash-dotted curve indicates the boundary of theblue region in the limit m / Λ → ∞ . The exclusion regions change only weakly with θ .The blue region is affected by the value of m / Λ while the screening constraint (green)remains the same. The self-accelerating cosmological models lie along the dashed lines. Theslope of these are always the same ( β = − α + const . ) but the intersection with the verticalaxis moves downwards if m / (2Ω Λ ) increases or if θ increases.in the infinite future, and α ∈ R and β ∈ R determine (among other things) how and if thescreening mechanism is active. The mixing angle can be expressed in terms of the ratio ofthe gravitational constants and the conformal factor,tan θ ≡ κc . (3.4)In the limit θ →
0, the physical metric coincides with the massless spin-2 field and we recovergeneral relativity (GR). In the limit θ → π/
2, the physical metric coincides with the massivespin-2 field and we recover dRGT (de Rham–Gabadadze–Tolley) massive gravity [51, 52]with a fixed auxiliary metric, see for example [53, 54].
Here, we summarize the analytical constraints on the physical parameters presented inRef. [38] and summarized in Fig. 1. To have a working screening mechanism restoring GRresults on local length scales (e.g., in the solar system), we must impose constraints on α and β , the exact form depending on θ and m FP . To have a continuous, real-valued back-ground cosmology devoid of the Higuchi ghost, we must put additional constraints on α and β , the exact form depending on θ and m / Λ . It is possible to express these constraintsin analytical form [38]. In particular, all B -parameters (except B and B ) must be nonzero.– 4 – .2 Special models In general, the physical parameters are independent. However, it is common to considersubsets of models where one or several of the B -parameters are set to zero, in which casethe physical parameters are not independent [39]. To be consistent with the analytical con-straints, B , B , and B must all be non-vanishing, excluding all two-parameters models (i.e.,with all but two of the B -parameters set to zero). The possible submodels are: B B B B (self-accelerating), B B B B , and B B B (minimal model). Here, we do not discuss the B B B B model explicitly, for two reasons: (i) it does not give as good fit to data as theself-accelerating model which has the same number of free parameters (see Section 5) and(ii) it exhibits the cosmological constant B which makes it less interesting from a theoreti-cal perspective. To connect with earlier studies, although incompatible with the analyticalconstraints, we study some two-parameter models in Appendix C. Self-accelerating models ( B B B B ). Cases where B = 0 are referred to as self-accelerating models and are of particular interest since there appears no cosmological constantterm in the general equations of motion for the physical metric g µν , see for example [38]. Still,in the equations of motions for the background cosmology there appears an effective cosmo-logical constant which is due to the interaction between the massless and massive spin-2fields. Setting B = 0 in (3.3), we get a linear relation between α and β , α + β − θ Ω Λ m , B = 0 , (3.5)that is, self-accelerating cosmologies lie along the line in the αβ -plane shown in Fig. 1.In order to have viable self-accelerating solutions satisfying the constraints introduced inSection 2, θ is bounded from above. An approximate upper limit is obtained by requiringthat the line (3.5) lies above the point ( α, β ) = ( − / , θ (cid:46) Λ /m , B = 0 , (3.6)which is an upper bound on θ for fixed values of m / (2Ω Λ ). Together with the Higuchibound m > Λ , this leaves us with a viable region close to m FP ∼ θ is very small,see Fig. 2. Thus, due to the Higuchi bound, there is no dRGT massive gravity limit ( θ → π/ m FP → ∞ enforces the GR limit θ → Λ ). Hence, for the self-acceleratingmodels to have interesting cosmological background solutions (i.e., different from ΛCDM), m FP cannot take too large values. From Fig. 1, we see that the limits α → ∞ and β → ∞ take us out from the allowed region and are not consistent with self-accelerating cosmologies. B B B (minimal) models. The B B B models are self-accelerating models where thecosmological constant terms in the two metric sectors are set to zero, that is B = B = 0,where B is the cosmological constant of the f µν metric. This is the most minimal modelcompatible with a working screening mechanism and a real-valued, continuous cosmologywithout the Higuchi ghost. For these models, a small value of the effective dark energydensity is technically natural in the sense of ’t Hooft [9]. Since B = B = 0, we can use The B B model (i.e., with only B and B non-vanishing) is still allowed. However, this is just twoindependent copies of general relativity. – 5 – igure 2 : Exclusion plot due to the analytical constraints (i.e., requiring a working screeningmechanism and a real-valued cosmology devoid of the Higuchi ghost). Here, we have setΩ Λ = 0 .
7. The green constraint applies to self-accelerating (i.e., B B B B ) models and theblue constraint to the B B B (minimal) models. For these, the viable region is a limitedarea close to m FP ∼
1. Remember that m FP is measured in units of H .(3.3) to express α and β in terms of θ and m / Ω Λ , α = − m / Ω Λ + 3 cot θ − θ m / Ω Λ , B = B = 0 (3.7a) β = 3 − m / Ω Λ + csc θ + 3 sec θ m / Ω Λ , B = B = 0 . (3.7b)A working Vainshtein mechanism requires β (cid:38) θ (cid:46) (cid:18) − m / Ω Λ ± (cid:113) m / Ω Λ ( − m / Λ ) (cid:19) , (3.8)see Fig. 2. Here, we must choose the minus sign corresponding to a positive α . The plussign implies a value for α not allowed by the analytical constraints. The largest value of θ isobtained by setting m = 2Ω Λ (from the Higuchi bound, see [38]). The result is, θ (cid:46) ◦ . (3.9)For the mixing angle not be close to the GR limit, the ratio m / Λ must be close to theHiguchi bound. This is indeed the case, see Tab. 1.– 6 – Local tests (stars, galaxies, and galaxy clusters)
Static, spherically symmetric (SSS) solutions can be used to approximate the gravitationalpotentials of for example the solar system, galaxies, and galaxy clusters. We refer to themas local solutions as opposed to cosmological solutions. In Refs. [12, 46], constraints wereplaced on the theory parameters assuming α ∼ β ∼
1. In Ref. [13], the authors constructeda phenomenological model, parameterized by the mixing angle θ and the graviton mass m FP .Here, we consider a general bimetric model including α and β . In principle, the effectivecosmological constant Ω Λ dictates the asymptotic structure at radial infinity. However, localsolutions are not sensitive to the value of Ω Λ which becomes influential only when approachingthe Hubble length scale, so Ω Λ does not affect the results of this section.There are two types of approximate analytical solutions to the SSS equations of mo-tion, applicable under different circumstances: the linearized solutions and the nonlinearVainshtein screening solutions. A conservative estimate is that the linearized solutions arecompatible with solar system data only if m FP (cid:38) or θ (cid:46) − [46]. In both cases, thebackground cosmology reduces to ΛCDM for all practical purposes and hence we loose manyof the interesting features of bimetric gravity [38]. Therefore, to ensure that our modelsare compatible with local tests of gravity while at the same time having novel cosmologi-cal solutions, we demand the existence of a working Vainshtein mechanism. It should bestressed however, that a careful analysis may reveal regions in the parameter space in whichthe linearized solutions pass the local tests even if m FP < or θ > − .The Vainshtein screening mechanism restores general relativity inside some radius, pro-vided that we satisfy the analytical constraints of Section 3.1; α and β determine howthe screening mechanism is realized (together with θ , m FP , and the mass of the source)[11, 12, 38, 55]. The radius within which GR can be restored is proportional to the Vain-shtein radius r V , r V ≡ (cid:0) M/H m (cid:1) / , (4.1)where M is the mass of the source. In the regions far inside and far outside the Vainshteinradius, the gravitational potential takes a simple form,Φ = − Mr , r (cid:28) r V , (4.2a)Φ = − Mr (cid:18) θ (cid:19) , r (cid:29) r V , (4.2b)and in the intermediate regions, it is calculated according to Ref. [38]. If α ∼ β ∼
1, theVainshtein radius sets the length scale within which we start to approach GR. As an example,to have a viable background cosmology, typically m FP ∼
1, and hence, r V /r ∗ ∼ ( ρ ∗ /ρ c ) / , (4.3)where r ∗ is the radius of the source and ρ ∗ is the mean density. Hence, the Vainshteinmechanism is relevant for all objects with a density greater than the critical density ρ c ≡ H /κ g and of course also all astronomical objects like the solar system, galaxies, and galaxy To be precise, the background cosmology of the large graviton mass limit reduces to ΛCDM, providedthat we do not impose any independent observation of the dimensionless matter density Ω m, . If such a valueis imposed, we also need to choose θ → – 7 – igure 3 : Exclusion plot, assuming α ∼ β ∼
1. Order of magnitude observational constraintsfrom galactic tests (orange) [12], solar system tests (purple) [46], galaxy cluster lensing (blue)[13], and gravitational waves (green) [45]. The Higuchi bound is shown in gray (here, we haveset Ω Λ = 0 .
7) [56]. Remember that m FP is measured in units of H .clusters, invalidating the applicability of linear structure formation severely for bimetricgravity [26]. For the Sun, the Vainshtein radius is r V ∼ AU where the gravitational forceis negligible anyway. For the Milky way, r Milky Way V ∼ r ∗ where we used M ∼ M (cid:12) and r ∗ ∼
50 000 ly.If we push α or β away from unity, we increase the radius within which we start toapproach GR [38]. Hence, α → ∞ and β → ∞ are GR limits for the local Vainshtein screeningsolutions. This is true for general bimetric models where all the physical parameters areindependent although, as discussed, inconsistent for self-accelerating models. The cases θ → m FP → α and β , the constraints placed on θ and m FP , assuming α ∼ β ∼ α ∼ β ∼
1, together with the Higuchi bound. Comparing Figs. 2 and 3, wesee that self-accelerating models (including the B B B models) which satisfy the analyticalconstraints automatically pass the local tests of gravity. The observational constraints thatwe use in this section are order of magnitude estimates. Detailed results can be found in thecorresponding references. The critical density is the total density (excluding curvature Ω k ) that the Universe must have today inorder to be spatially flat/Euclidean. – 8 – igure 4 : Confidence contours in the αβ -plane. The effective cosmological constant does notappear in the local solutions. For large α or β , GR is restored. The gray region is excludedby the requirement of a working Vainshtein mechanism, cf. Fig. 1. Left panel : Solar systemtests with θ = 1 . ◦ and m FP = 10 . Right panel : Gravitational lensing by a 10 M (cid:12) galaxy.Here, ( θ = 84 ◦ , m FP = 10 ). A conservative estimate is that bimetric gravity is compatible with solar system tests (at1 AU) except if θ (cid:38) − in the range 10 (cid:46) m FP (cid:46) [46]. The possibility of alleviatingthe constraints on θ and m FP by large values of α or β is illustrated in Fig. 4. There, weplot confidence contours in the αβ -plane with θ = 1 . ◦ and m FP = 10 which would beobservationally excluded if α ∼ β ∼
1. As a rough estimate, solar system tests constrainthe gravitational force to be proportional to 1 /r to an accuracy of 10 − [57]. Hence, wecalculate the chi-squared values as, χ = (cid:18) (Φ GR − Φ) / Φ GR − (cid:19) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) r =1 AU . (4.4)Here, Φ GR = − M/r . The confidence contours are surface levels of ∆ χ ≡ χ − χ where the 90 % level is at ∆ χ = 4 .
61 and the 95 % level is at ∆ χ = 5 .
99. For the best-fit cosmological parameters (see Section 5), the relative difference between Φ and Φ GR is (cid:39) − , so it is well within the observational constraints. As a rough model of gravitational lensing on galactic scales, we fix the galaxy mass andcalculate the lensing radius using the equation (in geometrized units), r lens = 4 M D source /R, (4.5)where D source is the distance of the source (which is set to 1 Gpc). For a galaxy of mass M = 10 M (cid:12) , the lensing radius is r lens (cid:39) ϕ ≡ (Φ + Ψ) / γ = Φ /ϕ. (4.6)We conservatively assume that measurements constrain γ to be within 10 % of the GR value γ GR = 1 [58], calculating the chi-squared values as, χ . lens . = (cid:18) γ ( r lens ) − . (cid:19) . (4.7)Assuming α ∼ β ∼
1, this constrains θ (cid:46) ◦ in the graviton mass range 40 (cid:46) m FP (cid:46) × [12]. Increasing α and β , these constraints are alleviated, see Fig. 4. For the best-fitcosmological parameters (see Section 5), the deviation of the gravitational slip from unity isat the level 10 − which is well within the theoretical constraints. In Ref. [13], the authors considered a phenomenological model for the local solutions, pa-rameterized by the mixing angle and the graviton mass. Comparing with clustering dataconstrains θ (cid:46) ◦ in the graviton mass range 10 (cid:46) m FP (cid:46) . As in the case of galactictests, the constraint can be alleviated by pushing α or β to large values. Since we have two metrics, the geometry of a homogeneous and isotropic universe is describedby two scale factors: a ( t ) for the physical metric g µν and (cid:101) a ( t ) for the auxiliary metric f µν .In the asymptotic future, the metrics tend to a proportional de Sitter solution f µν = c g µν ( c = const). Instead of the scale factor (cid:101) a , it is convenient to use the ratio between the scalefactors (including c and omitting the argument t ), y ≡ (cid:101) a/ ( ca ) , (5.1)which can be expressed in terms of the matter density Ω m by solving the quartic polynomial, −
13 cos θ m (1 + 2 α + β ) + (cid:20) Ω m + Ω Λ + m (cid:18) cos θ ( α + β ) − sin θ (cid:18) α + β (cid:19)(cid:19)(cid:21) y + m (cid:2) − cos θ β + sin θ (1 + 2 α + β ) (cid:3) y − (cid:20) Ω Λ + 13 m (cid:0) cos θ ( − α − β ) + 3 sin θ ( α + β ) (cid:1)(cid:21) y + 13 sin θ m βy = 0 . (5.2)The modified Friedmann equation for the scale factor a ( t ) reads (omitting the argument t ), E = Ω m + Ω k + Ω DE , E ≡ H/H . (5.3)Here, H = ˙ a/a is the Hubble parameter of the physical metric and H is the Hubble parameterevaluated today (i.e., the Hubble constant). Ω m and Ω k are dimensionless contributions frommatter and spatial curvature, respectively,Ω m = κ g ρ m H = Ω m, (1 + z ) w m ) , Ω k = − kH a = Ω k, (1 + z ) , (5.4)– 10 – .01 0.10 1 10 100 1000 - z Ω D E z R a t i oo f s c a l e f a c t o r s y Figure 5 : Examples of Ω DE ( z ) and y ( z ) for different physical parameters, Θ =( θ, m FP , Ω Λ , α, β ). Blue curve: Θ (cid:39) (18 ◦ , . , . , , . (cid:39) (18 ◦ , . , . , − . , . (cid:39) (0 . ◦ , , . , , (cid:39) (45 ◦ , , . , , DE and y decrease monotonically with redshift. The blue curve isa self-accelerating model, hence the dark energy density is zero in the early universe (i.e., as z → ∞ ). The red and green curves are models with positive dark energy density as z → ∞ while the purple model has negative dark energy density in the z → ∞ limit. The late-timecosmological constant is set by Ω Λ . In the limit z → ∞ , the dark energy density for thepurple model approaches B / (cid:39) − ρ m being the physical matter energy density. For simplicity, we consider only mattercontent with the same equation of state w m in eq. (5.3). However, we can generalize straight-forwardly by adding new terms for additional matter fields. Here, we are mostly interestedin redshifts 0 ≤ z (cid:46) w m = 0 (pressureless dust) in this range. In thelimit z → ∞ , we set w m = 1 / DE is a dynamical “dark energy” contributiondue to the gravitational interaction of the massive spin-2 field and can be expressed in termsof y as, Ω DE = Ω Λ − sin θ m (1 − y ) (cid:20) α (1 − y ) + β − y ) (cid:21) . (5.5)Since the equation for y (5.2) is a quartic polynomial, it has a closed-form solution with up tofour real solutions. However, only one of them is well-behaved. This is the finite/expandingbranch solution which is defined as the lowest lying, strictly positive root of y [39]. For thesesolutions, y and Ω DE increase monotonically with time (i.e., decrease with redshift), with y starting at y | z = ∞ = 0 at the Big Bang and ending up at y | z = − = 1 in the infinite future. Thedark energy Ω DE starts at some value at the Big Bang (negative, positive, or zero dependingon the value of the physical parameters) and ends up at the (positive) value Ω Λ in the infinitefuture, see Fig. 5 for some examples. For a general bimetric model, the equation of statefor the dark energy, w DE , starts at w DE | z = ∞ = − w DE | z = − = − w DE | z = ∞ = − (2+ w m ).)Hence, the bimetric fluid acts as a cosmological constant both in the early and late universe.However, since Ω DE is increasing with time, in the late universe it has a larger value than inthe early universe, see for example Fig. 5. In the intermediate region between the early-timeand late-time cosmological constant phases, there is a dynamical phase where the massivespin-2 field is dynamical, which can give rise to many different scenarios, depending on thevalues of the physical parameters, see Fig. 6 for some examples. Typically, in the intermediate– 11 – .01 0.10 1 10 100 1000 - - -
10 Redshift z w D E Θ =( ° ,10 ,0.7,10,10 ) Θ =( ° ,1.2,0.7,1,10.3 ) Θ =( ° ,2,0.75, - ) Θ =( ° ,2.3,0.74, - ) Figure 6 : Examples of the equation of state (EoS) for the bimetric dark energy fluid, fordifferent physical parameters Θ = ( θ, m FP , Ω Λ , α, β ). The red curve is a self-acceleratingmodel. All the models approach w DE = − w DE = −
1, except self-accelerating models which approach w DE = − (2 + w m ).In the intermediate region where the massive spin-2 field is dynamical, the equation of statedepends on the physical parameters. For the blue curve, Ω DE is negative in the early universe,hence it diverges when crossing zero at z (cid:39) .
3. Note that there is no physical singularityat this point. For the blue curve, the massive spin-2 field contributes as a dark mattercomponent (i.e., w DE = 0) in the redshift range 3 (cid:46) z (cid:46) w DE = − w DE = − (cid:46) z (cid:46)
10. For the green curve, the bimetric fluid contributes with a cosmological constantterm except in the range 1 (cid:46) z (cid:46)
100 where it has a phantom dark energy phase. The purplemodel is similar to the green one, with the difference that it also has a very brief phantomphase at z (cid:39) .
03. Note that equation of state does not diverge here but has a rather sharpdip reaching a minimum of w DE , min (cid:39) − w DE → − To explore the observational viability and to find the best fit parameters, we construct adiscrete grid in the parameter space Θ. For each point on the grid, we asses whether theanalytical constraints of Section 3.1 are satisfied. If not, the likelihood is set to zero. Thereby,we guarantee a working screening mechanism and that the cosmology is continuous, real-valued, and devoid of the Higuchi ghost. If the point is not excluded by these constraints,the likelihood is calculated by fitting the cosmological model to data from CMB, BAO, SNIa,and a background independent measurement of Ω m, . Since there is no established frameworkfor treating structure formation in bimetric theory, we will be conservative and combine CMBand BAO data in a way which effectively cancels the dependence on the cosmology before z ∗ (cid:39) z ∗ (cid:39) z ∈ [0 . , . (cid:96) A from Planck 2018 [64] as calculated in [65]. For the type– 12 – MB / BAO + XCLSNIa + XCLAll0.55 0.60 0.65 0.70 0.75 0.800.00.20.40.60.81.0 Ω Λ ℒ / ℒ m a x Flat Λ CDM
Figure 7 : Normalized likelihood distribution for spatially flat ΛCDM models as a functionof the cosmological constant Ω Λ . The best-fit point for CMB/BAO+XCL is Ω Λ = 0 . Λ = 0 . Λ = 0 .
695 for which χ = 52 . z ∈ [0 . , . m, = 0 . ± . . (5.6)As shown in Ref. [38], in the large parameter limits m FP , α, β → ∞ , the current matterdensity is given by, Ω m, = (1 − Ω Λ ) / cos θ, m FP , α, β → ∞ . (5.7)Hence, assuming Ω Λ (cid:39) .
7, eq. (5.7) together with (5.6) set an upper limit on θ in these largeparameter limits. Observations support Ω k (cid:39)
0, not only for a ΛCDM model but also forbimetric models [41]. Therefore, we assume a spatially flat universe with Ω k = 0, leaving uswith the bimetric physical parameters and Ω m, as cosmological parameters. (Since we onlyneed the expansion history up to z ∗ (cid:39) w m = 0.)Evaluating the modified Friedmann equation (5.3) and the y polynomial (5.2) today, we gettwo equations, containing Ω m, , y and the physical parameters. We solve the equations for y and Ω m, , leaving only the bimetric physical parameters independent. Hence, when fittingto cosmological data we are directly exploring the likelihood of the physical parameter space. Λ CDM.
In the spatially flat ΛCDM models with current matter density Ω m, = 1 − Ω Λ ,the cosmological constant Ω Λ is the only independent parameter. The result of fitting toCMB/BAO, SNIa, and Ω m, (from XCL) data is presented in Fig. 7. The value is Ω Λ =0 . ± .
014 (corresponding to Ω m, = 0 . ± . χ = 52 . B B (minimal) models. These models are a subset of the self-accelerating modelsand are of special interest since the cosmological constant terms B and B are absent. Still,the bimetric energy density Ω DE is dynamical and contributes to the accelerated expansion ofthe Universe. In Fig. 8a we show the two-dimensional marginalized confidence contours andparameter likelihoods of the B B B models. The two-dimensional marginalized confidencecontours are defined as level curves of ∆ χ ≡ χ − χ where 90 % corresponds to ∆ χ = 4 . χ = 5 .
99. To constrain the values of the physical parameters, we define theone-dimensional 90 % confidence interval as that within which ∆ χ < .
71. In Fig. 9 wecompile the marginalized 95 % confidence contours in the m FP θ -plane with the analyticalconstraints for general models, self-accelerating models, and B B B models.The massive spin-2 field effectively provides a dynamical phantom dark energy com-ponent which is influential in the range 0 (cid:46) z (cid:46)
10, see Fig. 10. The minimal χ value is χ = 52 .
2, hence there is a slight improvement in the fit compared with flat ΛCDM. Theconstraints on the physical parameters from these data sets and the analytical constraintsare shown in Tab. 1.The mixing angle θ has an upper bound θ < ◦ (90 % confidence) which is consistentwith the approximate analytical bound θ (cid:46) ◦ (3.9), see also Fig. 2. Remember that themixing angle parameterizes the proportions of the massive and massless gravitons in thephysical metric. Interestingly, there are models with a substantial mixing angle that stillhave viable background cosmology and local solutions (remember that our parameter choiceensures a working screening mechanism). The graviton mass has a lower limit m FP (cid:38) . m > Λ plus the fact that theeffective cosmological constant has a lower bound, Ω Λ = 0 . +0 . − . (90 % confidence).There is a degeneracy between θ and Ω Λ , see the “banana-like” shape in the Ω Λ θ -planeof Fig. 8a. The dark energy density increases with time and Ω Λ is the cosmological constantin the final de Sitter phase. For large values of the mixing angle θ , we are far away fromthat phase today (i.e., at z = 0) and hence the cosmological constant in the de Sitter phasei greater than today. The current matter density follows a distribution centered aroundΩ m, = 0 .
31, see Fig. 8d. The best-fit parameter values and constraints are summarized inTab. 1 and the expansion history for the best-fit model is plotted in Fig. 10. In terms ofthe B -parameters, the best-fit point is ( B , B , B , B , B ) best fit (cid:39) (0 , . , − . , . ,
0) and O ( B best fit n ) = 1. Self-accelerating models ( B B B B ). Similar results hold for the general self-acceleratingmodels. A novel feature that appears however, is a brief phase around z (cid:39) .
03 where theequation of state has a sharp dip, see Fig. 10. Since we have one additional free parame-ter compared to the B B B models, the marginalized confidence contours are broadened,compare Figs. 8a and 8b and the constraints on the physical parameters are weaker, seeTab. 1. In terms of the B -parameters, the best-fit point is ( B , B , B , B , B ) best fit (cid:39) (0 , . , − . , . , − . O ( B best fit n ) = 1. For the expansion history of the best-fitmodel, see Fig. 10. For the two-dimensional confidence contours in the parameters θ , m FP ,and Ω Λ , see Fig. 8b. In Fig. 11 of Appendix A, we plot the confidence contours also in β . The combined constraints from the analytical constraints and cosmological data sets aresummarized in Tab. 1. – 14 – a) B B B (minimal) models. (b) Self-accelerating models ( B B B B ).(c) General models ( B B B B B ). B B B Self - accGeneral Λ CDM0.30 0.32 0.340.00.51.0 Ω m ,0 ℒ / ℒ m a x (d) Matter density today, all models. Figure 8 : (a)-(c): Two-dimensional (marginalized) confidence contours in the physical pa-rameters θ (mixing angle), m FP (graviton mass), and Ω Λ (effective cosmological constant)when fitting to data from CMB/BAO+SNIa+XCL and imposing the analytical constraintsof Section 3.1. The one-dimensional plots (on the diagonals) are normalized likelihoods, L / L max . General models ( B B B B B ). For the most general bimetric models, the marginalizedconfidence contours are broadened compared to the self-accelerating models (compare Figs. 8band 8c) and the constraints on the physical parameters are weaker (see Tab. 1). In terms of the B -parameters, the best-fit point is ( B , B , B , B , B ) best fit (cid:39) ( − . , . , − . , . , − . O ( B best fit n ) = 10. The expansion history for the best-fit model is plotted in Fig. 10and the two-dimensional confidence contours of θ , m FP , and Ω Λ are plotted in Fig. 8c.– 15 – m FP Ω Λ α β Ω m, χ DoFGeneral model 11 ◦ +17 ◦ − ◦ . + ∞− . . +0 . − . − . + ∞−∞ + ∞− . +0 . − . ◦ +21 ◦ − ◦ . + ∞− . . +0 . − . – 29 + ∞− . +0 . − . B B B ◦ +15 ◦ − ◦ . + ∞− . . +0 . − . – – 0 . +0 . − . . +0 . − . – – 0 . +0 . − . Table 1 : Constraints on the physical parameters (90 % confidence) from the analyticalconstraints of Section 3.1, combined with CMB/BAO+SNIa+XCL data. The graviton mass m FP is expressed in units of H = 100 h km / s / Mpc = 2 . h × − eV /c . DoF (degreesof freedom) denotes the number of free parameters. (The upper constraint on Ω Λ for thegeneral model is calculated by linear extrapolation of the one-dimensional likelihood in therange Ω Λ ∈ [0 . , . Figure 9 : Marginalized 95 % confidence contours (below the curves) for B B B , self-accelerating ( B B B B ), and general ( B B B B B ) bimetric models. The colored regionsare excluded due to the analytical constraints. Gray region: Higuchi bound, applies to allmodels. Green region: excluded region for self-accelerating solutions, as approximated by(3.6). Blue region: excluded region for B B B models as approximated by (3.8).The confidence contours in the parameters α and β are plotted in Fig. 12 in Appendix A.The combined constraints on the physical parameters from the analytical constraints andcosmological data sets are summarized in Tab. 1. In Fig. 9 we compile the marginalized 95 %confidence contours in the m FP θ -plane with the analytical constraints for general models,self-accelerating models, and B B B models.Interestingly, the best-fit χ value does not improve when going from the self-acceleratingmodels to the general bimetric models, see Tab. 1. With respect to the data sets employedhere, one may conclude that the self-accelerating models are preferred compared to the gen-eral bimetric models due to the smaller number of free parameters. In other words, the– 16 – est - fit models ( E B R - E G R ) / E G R GeneralSelf - acc B B B %- % Redshift z Best - fit models Ω D E GeneralSelf - acc B B B - - Redshift z GeneralSelf - acc B B B Redshift z R a t i oo f s c a l e f a c t o r s y Best - fit models Best - fit models w D E GeneralSelf - acc B B B - - - Redshift z Figure 10 : Best-fit models.
Upper left panel: comparing the best-fit bimetric expansionrate, E BR , with the best-fit ΛCDM expansion rate, E GR . The expansion best-fit modelsfollow the ΛCDM model closely, except in the range 0 . (cid:46) z (cid:46) Upper right panel :the dark energy density Ω DE . For the general and self-accelerating models, there is a steepincrease in Ω DE around z = 0 .
03 which is due to the dip in the dark energy equation ofstate at that redshift.
Lower left panel : ratio of the scale factors y . Lower right panel : Darkenergy equation of state w DE . There is a divergence in w DE at z (cid:39) DE crossing zero. There is no physical singularity at this point.best-fit self-accelerating model is, to good approximation, the best-fit cosmological model.The results we have presented are sensitive the data sets employed. For example, if we onlyconsider statistical errors in the SNIa data, there is a more pronounced peak in the likelihoodof the physical parameters and the the best-fit is far away from any GR limit. Increasing theratio of the sound horizons at the drag epoch ( z d (cid:39) z ∗ (cid:39) R at this level, the peak in the likelihood is still far awayfrom any GR limit but is less pronounced than without R . To summarize, the observationalconstraints on the physical parameters depend on what data sets that we use. Here, wehave followed a very conservative approach and used the CMB/BAO ratio, excluding the R parameter (see Appendix B for more details).– 17 –CDM (flat) General bimetric Self-accelerating B B B AIC 54.5 62.0 60.0 58.2BIC 56.4 71.7 67.7 64.0∆AIC 0 +7 . . . . . . Table 2 : Information criteria. Here, we have set ΛCDM as our reference model so that∆IC = IC model − IC ΛCDM . A ∆AIC > +5 is commonly regarded as strong preference for thereference (flat ΛCDM) model compared to the (bimetric) model [69]. Information criteria.
As evident from Tab. 1, bimetric cosmology and ΛCDM provide,to close approximation, equally good fits to our data sets whereas the former has a largernumber of parameters to accomplish this. As a rough estimate of how well bimetric cosmologyperforms compared to ΛCDM, we compute the Akaike information criterion (AIC),AIC ≡ N param + χ , (5.8)and Bayesian information criterion (BIC),BIC = N param ln N data + χ , (5.9)where N param is the number of parameters in the model and N data is the number of datapoints [69]. The results are presented in Tab. 2. With a larger number of parameters,the performance, as assessed by the information criteria, is worse for the bimetric models,since flat ΛCDM provides such a good fit to these data sets. However, the theory mayhave beneficial properties that are not quantified by the information criteria, such as self-accelerating solutions.To summarize, bimetric cosmology and flat ΛCDM both provide good fits to our cos-mological data sets, with similar χ . Hence, the bimetric physical parameters are uncon-strained, except the effective cosmological constant and the upper bound on θ . At this level,the one-dimensional likelihoods of m FP , α , and β are approximately flat and the only con-straints on them are the analytical ones, see Section 3.1. Thanks to a working screeningmechanism, the best-fit cosmological parameters are compatible with observations from solarsystem tests, strong lensing by galaxies, and galaxy cluster lensing (see the previous subsec-tions). Constraints from gravitational waves, which constrain θ (cid:46) ◦ in the graviton massrange 10 (cid:46) m FP (cid:46) (see Fig. 3), are also satisfied. Due to the presence of four newparameters, bimetric cosmology is a more flexible model than ΛCDM. Interestingly, sincethere are viable bimetric models also away from the GR limits, there is a possibility thatother data sets can break the draw between the two. In particular, taking into account localmeasurements of the Hubble constant introduces a discrepancy in the ΛCDM model. Witha more flexible background evolution, bimetric cosmology might be able to alleviate the ten-sion. Also, with a framework for analyzing cosmological perturbations/structure formation,it may be possible to put more stringent constraints on the physical parameters.– 18 – .4 The Hubble tension Let us assess very roughly how well we expect the best-fit bimetric models to perform inrelation to the Hubble tension. Since the expansion is close to the best-fit ΛCDM model forredshifts z (cid:38)
10 (see Fig. 10), we assume here that the sound speed at photon decouplingtakes the same value as for a ΛCDM model. Because the expansion rate of the best-fitbimetric model is smaller in the range 0 . (cid:46) z (cid:46)
10, the integral I ∗ = (cid:82) z ∗ dz/E ( z ) is greaterthan the ΛCDM case. Therefore, in order for (cid:96) A to stay the same, H must be greater,see (B.2). Calculating the integral for the best-fit bimetric and ΛCDM models, we get H = 1 . H ΛCDM0 for the general model and H = 1 . H ΛCDM0 for the self-acceleratingmodel. Thus, instead of the ΛCDM value H ΛCDM0 = 67 . / s / Mpc [70] we get a valuewhich is greater by 0 . H = 68 . / s / Mpc, and a value that isgreater by 0 . H = 68 . / s / Mpc, for the self-accelerating model, which is somewhatcloser to the value from the late-time probes (cid:39)
73 km / s / Mpc [71]. For the best-fit B B B model, the Hubble constant is to good approximation unchanged compared with the ΛCDMvalue. Due to the presence of four new parameters, bimetric cosmology is a more flexiblemodel than ΛCDM and it is possible that there are viable models which increase the valueof H significantly. Ref. [10] concluded that bimetric two-parameter models can increase theHubble constant somewhat but not completely alleviate the tension. However, such modelsare restrictive and do not comply with the constraints from demanding a working screeningmechanism and a consistent cosmology. In Ref. [41], the Hubble tension is studied for generalbimetric models. However, they use the shift parameter R as a data point, which introducesa number of assumptions, see Appendix B for a detailed discussion. Whether or not thetension can be alleviated is a question for future study. We have studied the observational constraints on the physical parameters of bimetric gravityusing CMB/BAO, SNIa, and a background independent measurement of Ω m, (XCL) tofit for the background cosmology. Using the analytical constraints of [38], we ensure a real-valued background cosmology devoid of the Higuchi ghost and a working Vainshtein screeningmechanism for the local solutions. These constraints can be expressed analytically in termsof the physical parameters. Thereby, we ensure that gravitational tests ranging from solarsystem scales to galaxy cluster scales, are satisfied. (Constraints from gravitational waveobservation are also compatible with our results.) The results are summarized in Figs. 8 and9 and Tab. 1. For a general bimetric model, there is an upper limit on the mixing angle, θ < ◦ (90 % credibility). Hence, the physical metric can contain a significant amount of themassive spin-2 field. There is a lower limit on the graviton mass, m FP > . /c , m FP > . h × − eV /c corresponding to a Compton wavelength of thesize of the observable universe. The lower bound is many orders of magnitude smaller thanthe mass of the standard model particles. The effective cosmological constant, that is thedark energy density in the asymptotic future, is Ω Λ = 0 . +0 . − . (90 % credibility) and thecurrent matter density Ω m, = 0 . +0 . − . (90 % credibility). Interestingly, the general bimetricmodels (five free parameters) and the self-accelerating models (four free parameters) provideequally good fits, suggesting that the latter models are preferred over the former by thesedata sets. In the self-accelerating models, the accelerated expansion is due to the dynamicalmassive spin-2 field, without a cosmological constant– 19 –ince the bimetric models improve the χ only somewhat, compared to the flat ΛCDMmodel, the likelihoods of the physical parameters m FP , α , and β are approximately flat andhence there are no bounds on these parameters, except the analytical constraints. However,the results depend on the data sets employed. For example, if we set a greater value ofthe ratio of the sound horizons at the drag epoch and photon decoupling ( z d (cid:39) z ∗ (cid:39) θ (cid:46) ◦ )opens the possibility to break the draw between bimetric cosmology and ΛCDM by introduc-ing new data sets. In particular, taking late-time measurements of the Hubble constant intoaccount, there appears a discrepancy in the ΛCDM model. A back-of-the-envelope calcula-tion indicates that the best-fit bimetric model increases the value of the Hubble constant by (cid:39) . Acknowledgments
Thanks to Angelo Caravano and Marvin L¨uben for many interesting discussions on thesubject and to Fracesco Torsello for comments on the manuscript. This research utilizedthe Sunrise HPC facility supported by the Technical Division at the Department of Physics,Stockholm University. Special thanks to Mikica Kocic for technical support.– 20 – igure 11 : Confidence contours and normalized likelihoods for the self-accelerating models.As expected, the likelihoods of m FP and β are approximately flat above the thresholds setby the analytical constraints, m FP (cid:38) β (cid:38)
1. The former is due to the Higuchi bound, m > Λ , and the latter is due to the requirement of a working Vainshtein screeningmechanism. Note that m FP and β behave similarly. A Marginalized likelihoods
Here, we show the confidence contours and likelihoods in the parameters ( θ, m FP , Ω Λ , β ) forthe self-accelerating models in Fig. 11 and in the parameters ( θ, m FP , Ω Λ , α, β ) for the generalmodels in Fig. 12.In Fig. 13, we show the result of fitting B B B models to another data set. Here, we useCMB/BAO data with (cid:96) A = 301 . ± .
15 [72] and the same four BAO points as in Ref. [70].The SNIa data points are still from Pantheon but without systematic errors. The ratio of thesound horizons at the drag epoch and at photon decoupling is set to r s ( z d ) /r s ( z ∗ ) = 1 .
03. Asevident from Fig. 13, here the bimetric model improves the fit substantially compared with theflat ΛCDM model and the likelihood has a pronounced peak at ( θ, m FP , Ω Λ ) = (21 ◦ , . , . igure 12 : Confidence contours and normalized likelihoods for the general bimetric models.As expected, the likelihood of α is flat and the likelihoods of m FP and β are approximatelyflat above the thresholds set by the analytical constraints, m FP (cid:38) β (cid:38)
1. Thus, m FP , α , and β are poorly constrained by these cosmological data sets. There is a region in the αβ -plane which is excluded due to the constraints from having a working screening mechanism,cf. Fig. 1. B Data sets
Background independent Ω m, . The Chandra X-ray observatory measures the X-ray gasfraction in galaxy clusters and thereby obtains a measurement of the ratio of the baryonicmatter density and the total (pressureless) matter density. Combined with constraints on thephysical baryonic density from big bang nucleosynthesis (BBN) Ω b h = 0 . ± . h = 0 . ± .
08, the total matter density is Ω m, = 0 . ± .
04 [67]. Instead of the constrainton Ω b h from BBN, we use a constraint from extragalactic dispersion of fast radio bursts,which is independent of the background cosmology [68]. This gives a somewhat greater value– 22 – igure 13 : Confidence contours and normalized likelihoods for B B B models with adifferent data set than in the main text. Here, (cid:96) A is from [72] and the BAO points from [70].The SNIa data is from Pantheon but without systematic errors and the ratio of the soundhorizons is set to r s ( z d ) /r s ( z ∗ ) = 1 . b and hence the total matter density and the error also increases,Ω m, = 0 . ± . . (B.1)We refer to this data point as XCL (X-ray galaxy cluster). CMB data.
It is common to represent the information contained in the CMB power spec-trum by a few parameters, including the shift parameters (cid:96) A and R (see e.g. [73, 74]). Theformer is the angular scale of the sound horizon at photon decoupling, (cid:96) A , is given by (seee.g. [72]), (cid:96) A = πD A ( z ∗ ) /r s ( z ∗ ) , (B.2)where D A ( z ∗ ) is the comoving angular diameter distance to the last scattering surface z ∗ , D A ( z ) = I ( z ) /H , (B.3)with I ( z ) being, I ( z ) ≡ (cid:90) z dz (cid:48) E ( z (cid:48) ) , Ω k = 0 , (B.4)and r s ( z ∗ ) is the sound horizon at photon decoupling z ∗ . We set z ∗ (cid:39) R represents the angular scale of theHubble horizon at z ∗ [76], and is commonly given by, R = (cid:113) Ω m, H D A ( z ∗ ) . (B.5)– 23 –or a w CDM model, with the 2018 Planck data, the shift parameters are constrained to be[65], (cid:96) A = 301 . +0 . − . , R = 1 . +0 . − . , (B.6)with normalized covariance matrix, (cid:98) C = (cid:18) . .
47 1 (cid:19) . (B.7)The covariance matrix is given by C ij = σ i σ j (cid:98) C ij (no summation over the indices). In somepapers on bimetric cosmology, the CMB data is represented by (cid:96) A and R , see e.g. [10, 40, 41].To use R , one must introduce a number of assumptions. First, eq. (B.5) is applicable onlyif all contributions to E ( z ∗ ) are negligible, except the matter density Ω m ( z ∗ ) [77]. This iscertainly not the case if m FP → ∞ , α → ∞ , or β → ∞ . In these limits, the bimetricdark energy fluid mimics a cold dark matter component [38] and hence we have a “darkdegeneracy” [78]. In these cases, even if the background expansion history coincides with aΛCDM concordance model, R| m FP ,α,β →∞ = (cid:113) (1 − Ω Λ ) H D A ( z ∗ ) / cos θ. (B.8)Due to the cos θ in the denominator, R will be far away from the tabulated value even if thebackground expansion follows the best-fit ΛCDM model. Second, the shift parameters are notdirect observables, but are derived from the CMB power spectrum assuming a cosmologicalmodel (typically ΛCDM or w CDM), thus introducing assumptions of the primordial powerspectrum and the growth of structure [79]. In fact, R is more sensitive to these assumptionsthan (cid:96) A [79]. Since there is no established framework for treating structure formation inbimetric theory, we are conservative and combine (cid:96) A with BAO data in a way which effectivelyeliminates the dependence on the cosmology before z ∗ (cid:39) BAO data.
Baryon acoustic oscillations measures the ratio of the sound horizon scale r s ( z d ) at the drag epoch z d and a particular cosmological distance scale, here D A , d A , or D V (depending on the data set), collectively denoted by D X . D A is the comoving angulardiameter distance (B.3), d A is the angular diameter distance, d A ( z ) = D A ( z ) / (1 + z ) , (B.9)and D V is the volume average distance, D V ( z ) = (cid:2) D A ( z ) z/H ( z ) (cid:3) / . (B.10)We use the data sets from 6dFGS [59], SDSS MGS [60], BOSS DR12 [61], BOSS DR14 [62],and eBOSS QSO [63], in total ten points in the redshift range z ∈ [0 . , . i = π D X ( z i ) /r s ( z d ) (cid:96) A r s ( z d ) r s ( z ∗ ) = D X ( z i ) D A ( z ∗ ) , (B.11)where i runs over the BAO points z i ∈ { . , ..., . } and D X is a cosmological distancescale depending on the data set, according to Tab. 3. The drag epoch z d is the time where– 24 –ata set z eff Distance measure 10 (cid:98) C ij . D V /r d = 2 . ± .
133 —SDSS MGS [60] 0 . D V /r d = 4 . ± .
168 —0.38 D A /r d = 10 . ± .
15 10 . D A /r d = 13 . ± .
18 4970 10 . D A /r d = 15 . ± .
22 1991 984 10 BOSS DR14 [62] 0 . D V /r d = 16 . ± .
41 —0 . d A /r d = 10 . ± . . d A /r d = 12 . ± . . d A /r d = 11 . ± .
65 2662 6130 10 . d A /r d = 12 . ± .
99 248 954 4257 10 Table 3 : BAO data sets used in our analysis. Here, r d = r s ( z d ) and the covariance matrixis given by C ij = σ i σ j (cid:98) C ij (no summation implied).the baryons were released from the Compton drag of the photons and hence the acousticoscillations frozen in. From the equations in Ref. [81], we calculate z d (cid:39) z ∗ and thus Π i depends only on the expansion history between photon decouplingand today [80]. The ratio of the sound horizons at z d (cid:39) z ∗ (cid:39) r s ( z ∗ ) = (144 . ± .
48) Mpc and r s ( z d ) = (147 . ± .
48) Mpc [64], so that, r s ( z d ) /r s ( z ∗ ) = 1 . ± . . (B.12)We calculate the right-hand side of (B.11), Π model i , for each bimetric model and the valueof the left-hand side, Π obs i , is obtained through the observed quantities (B.6) and (B.12)together with the Tab. 3. The χ for each model is calculated as, χ / BAO = (cid:88) ij (Π model i − Π obs i ) C − ij (Π model j − Π obs j ) . (B.13) SNIa data.
Since the properly calibrated luminosity of supernovae of type Ia is believedto be independent of redshift, they can used as standard candles, allowing us to calibratecosmological distances. Here, we make use of binned data of the peak apparent magnitude m B from the Pantheon data set; 40 bins ranging from z = 0 .
014 to z = 1 .
61 [66]. m B isgiven by, m B = M + 5 log d L , (B.14) Available here, last checked 2020-11-05. – 25 –here d L is the luminosity distance rescaled by a factor of H /c (here, c is the speed of light), d L = (1 + z ) (cid:90) z dz (cid:48) E ( z (cid:48) ) , (B.15)and M is, M = 25 + M B + 5 log ( c/H ) , (B.16)with M B being the absolute magnitude of the supernova. Here, we marginalize over M , seefor example [82], with the result, χ = (cid:98) χ − B /C + ln( C/ π ) , (B.17)where, ∆ i ≡ d L ( z i ) − m B ( z i ) , (cid:98) χ ≡ (cid:88) ij ∆ i C − ij ∆ j , (B.18) B ≡ (cid:88) ij u i C − ij ∆ j , C ≡ (cid:88) ij u i C − ij u j , (B.19)and u i = 1 , i ∈ { , , ..., } . (B.20) C ij is the covariance matrix, including both statistical and systematic errors. C Two-parameter models
Here, all B -parameters but two are set to zero. To have a real-valued cosmology in the earlyuniverse, B > B and B . However, we include them here to connect toearlier works. As an example, in a B B model (i.e., B = B = B = 0), we can let θ andΩ Λ be our independent parameters. Then α , β , and m FP can be expressed in terms of theseas, m = 12 (csc θ + 3 sec θ )Ω Λ , B = B = B = 0 , (C.1a) α = 1 − θ θ , B = B = B = 0 , (C.1b) β = 0 , B = B = B = 0 . (C.1c)The B -parameters read, B = 32 (1 − tan θ )Ω Λ , B = B = B = 0 , (C.2a) B = −
12 (1 − θ )Ω Λ , B = B = B = 0 . (C.2b)For more examples, see [39]. A similar exercise can be carried out for the other self-accelerating two-parameter models B B and B B .Assuming that the best fit to these models is close to a ΛCDM model, we can predictanalytically the observationally viable regions in the parameter space, without making anyactual data fitting. Setting Ω Λ = 0 .
7, the GR limit θ → Λ (cid:39) .
7. Hence, we expect ourbest-fit models to lie in the neighborhood of this line, see Fig. 14.– 26 – igure 14 : Viable regions in the parameter spaces of self-accelerating two-parameter models: B B (top), B B (middle), and B B (bottom). The black curves are traced out by letting θ → Λ = 0 .
70. Right panel: results in the B -parameter plane, with B n definedas in (3.1). Left panel: the ¯ β -parameters can be identified with the B -parameters of [6] whichis the same as the β -parameters of [8]. Indeed, the analytical prediction is consistent with[6, 8] where the bimetric models are fitted to cosmological data. D Numerical details
D.1 Numerical errors
Numerical errors are controlled by estimating upper limits in the derived quantities. Ulti-mately, we are interested in the χ at each point of the parameter space, which we calculateusing (B.13) and (B.17), depending on E ( z ) and I ( z ) ≡ (cid:82) z dz (cid:48) /E ( z (cid:48) ). We estimate an upperlimit on the numerical errors of these. We introduce a discrete logarithmic grid (plus thepoint z = 0) containing 200 points in the relevant redshift range 0 ≤ z ≤ z ∗ .We start by evaluating the equations (5.2) and (5.3) at z = 0 (i.e., today) and solve for– 27 – and Ω m, . Rewriting the equations,Ω m, = 1 − Ω Λ + m sin θ (1 − y ) (cid:20) α (1 − y ) + β − y ) (cid:21) , (D.1a)0 = − Λ y + 13 cos θ m FP y (cid:2) α + β − α + β ) y + 3 βy − (1 − α + β ) y (cid:3) . (D.1b)From Ω m, , we know Ω m = Ω m, (1 + z ) w m ) as a function of redshift and y can be solvedat each redshift from eq. (5.2). Finally, E is computed by E ( z ) = (cid:112) Ω m ( z ) + Ω DE ( y ( z ))(5.3) at each z .When solving the equations for y , Ω m, , and y , we introduce numerical errors. Wesolve the equations to 15 digit precision, meaning that the relative errors in these quantitiesare 10 − . These numerical errors propagate to E via, (cid:20) ∆ EE (cid:21) numerical y = 12 y (cid:101) Ω DE ∂ (cid:101) Ω DE ∂y ∆ yy , (D.2a) (cid:20) ∆ EE (cid:21) numerical Ω m, = 12 Ω m (cid:101) Ω DE ∂ (cid:101) Ω DE ∂y (cid:18) ∂ Ω m ∂y (cid:19) − ∆Ω m, Ω m, , (D.2b) (cid:20) ∆ EE (cid:21) numerical y = 12 Ω m Ω m, y (cid:101) Ω DE ∂ (cid:101) Ω DE ∂y (cid:18) ∂ Ω m ∂y (cid:19) − ∂ Ω m, ∂y ∆ y y . (D.2c)Here, ∆ denotes the numerical error in a calculated quantity. We are now ready to show howthis error affect the χ value. Starting with CMB/BAO, the numerical error in χ comesfrom Π model i which has contributions from ∆ I ∗ , ∆ I ( z i ), and ∆ E (see eq. (B.13)). Using thechain rule and the triangle inequality, (cid:12)(cid:12)(cid:12)(cid:12) ∆Π i Π i (cid:12)(cid:12)(cid:12)(cid:12) ≤ (cid:12)(cid:12)(cid:12)(cid:12) ∆ I ∗ I ∗ (cid:12)(cid:12)(cid:12)(cid:12) + 23 (cid:12)(cid:12)(cid:12)(cid:12) ∆ I ( z i ) I ( z i ) (cid:12)(cid:12)(cid:12)(cid:12) + 13 (cid:12)(cid:12)(cid:12)(cid:12) ∆ EE (cid:12)(cid:12)(cid:12)(cid:12) , (D.3)where the last term is given in (D.2). In the integral terms, there are two types of errorscontributing: the numerical error in E and the finite redshift grid size. Concerning the firsttype, they tend to cancel when integrating. However, as an upper limit we estimate it to beof the same magnitude as ∆ E/E . The errors due to the finite redshift grid size is estimatedby comparing the value of the integral using a grid with 100 points and a grid with 200points. The difference in their value compared to the value of the integral gives an upperlimit on the numerical error, ∆
II < I
100 points − I
200 points I
200 points . (D.4)Using the variance formula, we finally get the upper limit on the numerical error for χ , (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ∆ χ / BAO χ / BAO (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) < χ / BAO
Max (cid:12)(cid:12)(cid:12)(cid:12) ∆Π i Π i (cid:12)(cid:12)(cid:12)(cid:12) (cid:88) i =1 (cid:12)(cid:12)(cid:12)(cid:12) Π i Π i − Π obs , i σ i (cid:12)(cid:12)(cid:12)(cid:12) (D.5)Similarly, for the SNIa calculations, (cid:12)(cid:12)(cid:12)(cid:12) ∆ χ χ (cid:12)(cid:12)(cid:12)(cid:12) < χ ln 10 (cid:88) i =1 σ i (cid:12)(cid:12)(cid:12)(cid:12) ∆ i − BC (cid:12)(cid:12)(cid:12)(cid:12) (cid:12)(cid:12)(cid:12)(cid:12) ∆ I ( z i ) I ( z i ) (cid:12)(cid:12)(cid:12)(cid:12) , (D.6)where i = 1 , ...,
40 and z , ..., z are the redshift of the different bins. The resulting upperlimits on the numerical errors are presented in Fig. 15.– 28 – elf - acc model B B B modelGeneral model - - - % % % Δχ / χ R e l a t i v e f r e q u e n c y Numerical errors ( upper bounds ) Figure 15 : The distribution of relative errors in χ (note that these errors are upper limits).In all cases, ∆ χ /χ < − . D.2 Scanning details
In the scanning process, we implement the constraints of Section 2 by assigning zero likelihoodto the points in parameter space which violate these constraints. A couple of numericalproblems can occur at each point in the parameter space. First, solving the equations for y and y , the numerical algorithm may not find the finite branch root, that is within the range0 < y <
1. Second, the algorithm may not be able to solve the equations to the requiredprecision (15 digits). Both these problems can in principle be solved by refining the numericalalgorithm. If any of these errors occur, we remove that point from the grid. To ensure thatthe final error in χ is small enough, we remove points with ∆ χ /χ > − .Scanning the general models (Fig. 8c), we use a grid of size 30 = 24 . × pointsof which (cid:39) . (cid:39) . × points whereof 0 . B B B models are a grid size of 80 (cid:39) . × points whereof none were removed. References [1] M. S. Volkov,
Cosmological solutions with massive gravitons in the bigravity theory , JHEP (2012) 035, [ ].[2] M. von Strauss, A. Schmidt-May, J. Enander, E. Mortsell and S. F. Hassan, CosmologicalSolutions in Bimetric Gravity and their Observational Tests , JCAP (2012) 042,[ ].[3] D. Comelli, M. Crisostomi, F. Nesti and L. Pilo,
FRW Cosmology in Ghost Free MassiveGravity , JHEP (2012) 067, [ ].[4] M. S. Volkov, Hairy black holes in the ghost-free bigravity theory , Phys. Rev.
D85 (2012)124043, [ ].[5] M. S. Volkov,
Exact self-accelerating cosmologies in the ghost-free massive gravity – the detailedderivation , Phys. Rev. D (2012) 104022, [ ].[6] Y. Akrami, T. S. Koivisto and M. Sandstad, Accelerated expansion from ghost-free bigravity: astatistical analysis with improved generality , JHEP (2013) 099, [ ]. – 29 –
7] M. S. Volkov,
Self-accelerating cosmologies and hairy black holes in ghost-free bigravity andmassive gravity , Class. Quant. Grav. (2013) 184009, [ ].[8] F. Koennig, A. Patil and L. Amendola, Viable cosmological solutions in massive bimetricgravity , JCAP (2014) 029, [ ].[9] G. ’t Hooft,
Naturalness, chiral symmetry, and spontaneous chiral symmetry breaking , NATOSci. Ser. B (1980) 135–157.[10] E. M¨ortsell and S. Dhawan, Does the Hubble constant tension call for new physics? , JCAP (2018) 025, [ ].[11] E. Babichev and M. Crisostomi,
Restoring general relativity in massive bigravity theory , Phys.Rev.
D88 (2013) 084002, [ ].[12] J. Enander and E. Mortsell,
On stars, galaxies and black holes in massive bigravity , JCAP (2015) 023, [ ].[13] M. Platscher, J. Smirnov, S. Meyer and M. Bartelmann,
Long Range Effects in GravityTheories with Vainshtein Screening , JCAP (2018) 009, [ ].[14] D. Comelli, M. Crisostomi and L. Pilo, Perturbations in Massive Gravity Cosmology , JHEP (2012) 085, [ ].[15] N. Khosravi, H. R. Sepangi and S. Shahidi, Massive cosmological scalar perturbations , Phys.Rev.
D86 (2012) 043517, [ ].[16] M. Berg, I. Buchberger, J. Enander, E. Mortsell and S. Sjors,
Growth Histories in BimetricMassive Gravity , JCAP (2012) 021, [ ].[17] Y. Sakakihara, J. Soda and T. Takahashi,
On Cosmic No-hair in Bimetric Gravity and theHiguchi Bound , PTEP (2013) 033E02, [ ].[18] F. K¨onnig and L. Amendola,
Instability in a minimal bimetric gravity model , Phys. Rev.
D90 (2014) 044030, [ ].[19] D. Comelli, M. Crisostomi and L. Pilo,
FRW Cosmological Perturbations in Massive Bigravity , Phys. Rev.
D90 (2014) 084003, [ ].[20] A. De Felice, A. E. G¨umr¨uk¸c¨uo˘glu, S. Mukohyama, N. Tanahashi and T. Tanaka,
Viablecosmology in bimetric theory , JCAP (2014) 037, [ ].[21] A. R. Solomon, Y. Akrami and T. S. Koivisto,
Linear growth of structure in massive bigravity , JCAP (2014) 066, [ ].[22] F. Koennig, Y. Akrami, L. Amendola, M. Motta and A. R. Solomon,
Stable and unstablecosmological models in bimetric massive gravity , Phys. Rev.
D90 (2014) 124014, [ ].[23] M. Lagos and P. G. Ferreira,
Cosmological perturbations in massive bigravity , JCAP (2014) 026, [ ].[24] F. K¨onnig,
Higuchi Ghosts and Gradient Instabilities in Bimetric Gravity , Phys. Rev.
D91 (2015) 104019, [ ].[25] K. Aoki, K.-i. Maeda and R. Namba,
Stability of the Early Universe in Bigravity Theory , Phys.Rev.
D92 (2015) 044054, [ ].[26] E. Mortsell and J. Enander,
Scalar instabilities in bimetric gravity: The Vainshtein mechanismand structure formation , JCAP (2015) 044, [ ].[27] Y. Akrami, S. F. Hassan, F. K¨onnig, A. Schmidt-May and A. R. Solomon,
Bimetric gravity iscosmologically viable , Phys. Lett.
B748 (2015) 37–44, [ ].[28] M. H¨og˚as, F. Torsello and E. M¨ortsell,
On the stability of bimetric structure formation , JCAP (2020) 046, [ ]. – 30 –
29] M. L¨uben, A. Schmidt-May and J. Smirnov,
Vainshtein Screening in Bimetric Cosmology , .[30] M. Kocic, Geometric mean of bimetric spacetimes , .[31] M. Kocic, Causal propagation of constraints in bimetric relativity in standard 3+1 form , JHEP (2019) 219, [ ].[32] M. Kocic, A. Lundkvist and F. Torsello, On the ratio of lapses in bimetric relativity , Class.Quant. Grav. (2019) 225013, [ ].[33] F. Torsello, M. Kocic, M. H¨og˚as and E. M¨ortsell, Covariant BSSN formulation in bimetricrelativity , Class. Quant. Grav. (2020) 025013, [ ].[34] M. Kocic, F. Torsello, M. H¨og˚as and E. Mortsell, Spherical dust collapse in bimetric relativity:Bimetric polytropes , .[35] F. Torsello, The mean gauges in bimetric relativity , Class. Quant. Grav. (2019) 235010,[ ].[36] F. Torsello, bimEX : A Mathematica package for exact computations in bimetric relativity , Comput. Phys. Commun. (2020) 106948, [ ].[37] M. Kocic, F. Torsello, M. H¨og˚as and E. M¨ortsell,
Initial data and first evolutions of dust cloudsin bimetric relativity , Class. Quant. Grav. (2020) 165010.[38] M. H¨og˚as and E. M¨ortsell, in preparation , .[39] M. L¨uben, A. Schmidt-May and J. Weller, Physical parameter space of bimetric theory andSN1a constraints , .[40] S. Dhawan, A. Goobar, E. M¨ortsell, R. Amanullah and U. Feindt, Narrowing down the possibleexplanations of cosmic acceleration with geometric probes , JCAP (2017) 040,[ ].[41] M. Lindner, K. Max, M. Platscher and J. Rezacek,
Probing alternative cosmologies through theinverse distance ladder , .[42] A. De Felice, T. Nakamura and T. Tanaka, Possible existence of viable models of bi-gravity withdetectable graviton oscillations by gravitational wave detectors , PTEP (2014) 043E01,[ ].[43] M. Fasiello and R. H. Ribeiro,
Mild bounds on bigravity from primordial gravitational waves , JCAP (2015) 027, [ ].[44] G. Cusin, R. Durrer, P. Guarato and M. Motta,
Inflationary perturbations in bimetric gravity , JCAP (2015) 043, [ ].[45] K. Max, M. Platscher and J. Smirnov, Gravitational Wave Oscillations in Bigravity , Phys. Rev.Lett. (2017) 111101, [ ].[46] M. L¨uben, E. M¨ortsell and A. Schmidt-May,
Bimetric cosmology is compatible with local testsof gravity , .[47] S. Sjors and E. Mortsell, Spherically Symmetric Solutions in Massive Gravity and Constraintsfrom Galaxies , JHEP (2013) 080, [ ].[48] J. Enander and E. M¨ortsell, Strong lensing constraints on bimetric massive gravity , JHEP (2013) 031, [ ].[49] D. G. Boulware and S. Deser, Can gravitation have a finite range? , Phys. Rev. D6 (1972)3368–3382.[50] S. F. Hassan and R. A. Rosen, Bimetric Gravity from Ghost-free Massive Gravity , JHEP (2012) 126, [ ]. – 31 –
51] C. de Rham and G. Gabadadze,
Generalization of the Fierz-Pauli Action , Phys. Rev.
D82 (2010) 044020, [ ].[52] C. de Rham, G. Gabadadze and A. J. Tolley,
Resummation of Massive Gravity , Phys. Rev.Lett. (2011) 231101, [ ].[53] S. F. Hassan, A. Schmidt-May and M. von Strauss,
On Consistent Theories of Massive Spin-2Fields Coupled to Gravity , JHEP (2013) 086, [ ].[54] S. F. Hassan, A. Schmidt-May and M. von Strauss, Particular Solutions in Bimetric Theoryand Their Implications , Int. J. Mod. Phys.
D23 (2014) 1443002, [ ].[55] A. I. Vainshtein,
To the problem of nonvanishing gravitation mass , Phys. Lett. (1972) .[56] A. Higuchi,
Forbidden Mass Range for Spin-2 Field Theory in De Sitter Space-time , Nucl.Phys. B (1987) 397–436.[57] C. M. Will,
The Confrontation between General Relativity and Experiment , Living Rev. Rel. (2014) 4, [ ].[58] J. Schwab, A. S. Bolton and S. A. Rappaport, Galaxy-Scale Strong Lensing Tests of Gravityand Geometric Cosmology: Constraints and Systematic Limitations , Astrophys. J. (2010)750–757, [ ].[59] F. Beutler, C. Blake, M. Colless, D. H. Jones, L. Staveley-Smith, L. Campbell et al.,
The 6dfgalaxy survey: baryon acoustic oscillations and the local hubble constant , Monthly Notices ofthe Royal Astronomical Society (Jul, 2011) 3017–3032.[60] A. J. Ross, L. Samushia, C. Howlett, W. J. Percival, A. Burden and M. Manera,
The clusteringof the SDSS DR7 main Galaxy sample –I. A 4 per cent distance measure at z = 0 . Mon.Not. Roy. Astron. Soc. (2015) 835–847, [ ].[61]
BOSS collaboration, S. Alam et al.,
The clustering of galaxies in the completed SDSS-IIIBaryon Oscillation Spectroscopic Survey: cosmological analysis of the DR12 galaxy sample , Mon. Not. Roy. Astron. Soc. (2017) 2617–2652, [ ].[62] J. E. Bautista et al.,
The SDSS-IV extended Baryon Oscillation Spectroscopic Survey: BaryonAcoustic Oscillations at redshift of 0.72 with the DR14 Luminous Red Galaxy Sample , Astrophys. J. (2018) 110, [ ].[63] G.-B. Zhao et al.,
The clustering of the SDSS-IV extended Baryon Oscillation SpectroscopicSurvey DR14 quasar sample: a tomographic measurement of cosmic structure growth andexpansion rate based on optimal redshift weights , Mon. Not. Roy. Astron. Soc. (2019)3497–3513, [ ].[64]
Planck collaboration, N. Aghanim et al.,
Planck 2018 results. VI. Cosmological parameters , Astron. Astrophys. (2020) A6, [ ].[65] L. Chen, Q.-G. Huang and K. Wang,
Distance Priors from Planck Final Release , JCAP (2019) 028, [ ].[66] D. Scolnic et al., The Complete Light-curve Sample of Spectroscopically Confirmed SNe Ia fromPan-STARRS1 and Cosmological Constraints from the Combined Pantheon Sample , Astrophys.J. (2018) 101, [ ].[67] S. Allen, R. Schmidt, H. Ebeling, A. Fabian and L. van Speybroeck,
Constraints on darkenergy from Chandra observations of the largest relaxed galaxy clusters , Mon. Not. Roy.Astron. Soc. (2004) 457, [ astro-ph/0405340 ].[68] J. Macquart, J. Prochaska, M. McQuinn et al.,
A census of baryons in the Universe fromlocalized fast radio bursts , Nature (2020) 391–395.[69] A. R. Liddle,
Information criteria for astrophysical model selection , Mon. Not. Roy. Astron.Soc. (2007) L74–L78, [ astro-ph/0701113 ]. – 32 – Planck collaboration, P. Ade et al.,
Planck 2015 results. XIII. Cosmological parameters , Astron. Astrophys. (2016) A13, [ ].[71] A. G. Riess et al.,
A 2.4% Determination of the Local Value of the Hubble Constant , Astrophys.J. (2016) 56, [ ].[72]
Planck collaboration, P. Ade et al.,
Planck 2015 results. XIV. Dark energy and modifiedgravity , Astron. Astrophys. (2016) A14, [ ].[73] A. Kosowsky, M. Milosavljevic and R. Jimenez,
Efficient cosmological parameter estimationfrom microwave background anisotropies , Phys. Rev. D (2002) 063007, [ astro-ph/0206014 ].[74] Y. Wang and P. Mukherjee, Observational Constraints on Dark Energy and Cosmic Curvature , Phys. Rev. D (2007) 103533, [ astro-ph/0703780 ].[75] W. Hu and N. Sugiyama, Small scale cosmological perturbations: An Analytic approach , Astrophys. J. (1996) 542–570, [ astro-ph/9510117 ].[76] G. Efstathiou and J. Bond,
Cosmic confusion: Degeneracies among cosmological parametersderived from measurements of microwave background anisotropies , Mon. Not. Roy. Astron. Soc. (1999) 75–97, [ astro-ph/9807103 ].[77]
WMAP collaboration, E. Komatsu et al.,
Five-Year Wilkinson Microwave Anisotropy Probe(WMAP) Observations: Cosmological Interpretation , Astrophys. J. Suppl. (2009) 330–376,[ ].[78] M. Kunz,
The dark degeneracy: On the number and nature of dark components , Phys. Rev. D (2009) 123001, [ astro-ph/0702615 ].[79] O. Elgaroy and T. Multamaki, On using the CMB shift parameter in tests of models of darkenergy , Astron. Astrophys. (2007) 65, [ astro-ph/0702343 ].[80] J. Sollerman et al.,
First-Year Sloan Digital Sky Survey-II (SDSS-II) Supernova Results:Constraints on Non-Standard Cosmological Models , Astrophys. J. (2009) 1374–1385,[ ].[81] D. J. Eisenstein and W. Hu,
Baryonic features in the matter transfer function , Astrophys. J. (1998) 605, [ astro-ph/9709112 ].[82] M. Goliath, R. Amanullah, P. Astier, A. Goobar and R. Pain,
Supernovae and the nature of thedark energy , Astron. Astrophys. (2001) 6–18, [ astro-ph/0104009 ].].