N-body Simulations of $\gamma$ Gravity

0.65 MB·

by  Marcelo Vargas dos Santos

Astronomy&Astrophysics587,A132(2016) (cid:13)cESO2016
March3,2016
N-body simulations of γ gravity
MarceloVargasdosSantos1,2,(cid:63),HansA.Winther3,(cid:63)(cid:63),DavidF.Mota2,(cid:63)(cid:63)(cid:63),andIoavWaga1,(cid:63)(cid:63)(cid:63)(cid:63)
1 InstitutodeFísica,UniversidadeFederaldoRiodeJaneiroC.P.68528,CEP21941-972,RiodeJaneiro,RJ,Brazil
2 InstituteofTheoreticalAstrophysics,UniversityofOslo,Postboks1029,0315Oslo,Norway
3 Astrophysics,UniversityofOxford,DWB,KebleRoad,Oxford,OX13RH,UK
March3,2016
6
1 ABSTRACT
0
Wehaveinvestigatedstructureformationintheγ gravity f(R)modelwithN-bodysimulations.Theγ gravitymodelisaproposal
2
which,unlikeotherviable f(R)models,notonlychangesthegravitationaldynamics,butcaninprinciplealsohavesignaturesatthe

r backgroundlevelthataredifferentfromthoseobtainedinΛCDM(Cosmologicalconstant,ColdDarkMatter).Theaimofthispaper
a istostudythenonlinearregimeofthemodelinthecasewhere,atlatetimes,thebackgrounddiffersfromΛCDM.Wequantifythe
M
signaturesproducedonthepowerspectrum,thehalomassfunction,andthedensityandvelocityprofiles.Toappreciatethefeatures
ofthemodel,wehavecomparedittoΛCDMandtheHu-Sawicki f(R)models.Fortheconsideredsetofparameterswefindthatthe
2 screeningmechanismisineffective,whichgivesrisetodeviationsinthehalomassfunctionthatdisagreewithobservations.Thisdoes
notruleoutthemodelperse,butrequireschoicesofparameterssuchthat|fR0|ismuchsmaller,whichwouldimplythatitscosmic
] expansionhistorycannotbedistinguishedfromΛCDMatthebackgroundlevel.
O
C Keywords. Gravitation–Cosmology:darkenergy–Galaxies:clusters:general–Cosmology:large-scalestructureofUniverse–
Galaxies:halos
.
h
p
-1. Introduction 1980),isdrivenbyatermofthetype f(R)=αR2 (α>0)andis
o
stillinexcellentaccordancewithobservations(Adeetal.2014).
r
Sincethediscoveryin1998(Riessetal.1998;Perlmutteretal.
t More recently, the idea of an acceleration driven by late-time
s1999)thattheUniverseisspeedingupinsteadofslowingdown
a curvaturehasalsobeenexploredinCapozzielloetal.(2003)and
(aswouldbeexpectedifgravityisalwaysattractive),consider-
[ Carrolletal.(2004).Theseauthorsconsideredatheoryinwhich
ableefforthasbeendevotedtounderstandingthephysicalmech- f(R)=−αR−n(n>0andα>0).However,thesemodelsdonot
2anism behind this cosmic acceleration. The two main theoreti-
have a regular matter-dominated era and are incompatible with
vcal approaches considered in the literature to explain this phe-
structureformation(Adeetal.2014).
3
nomena are (1) to assume the existence of a new component
3with a sufficiently negative pressure (p < −ρ/3), generically Tobuildacosmologicallyviable f(R)theory,somestability
4 conditions have to be satisfied (Pogosian & Silvestri 2008): (a)
5denoted dark energy, and (2) to consider that general relativity f ≡d2f/dR2 >0(notachyons);(b)1+f ≡1+df/dR>0[the
01.0hcdouaafsrtsevttaaoittsuebrEeewin(mmstoeod≡idnifii’fisepedcdoa/sgtρmrlaaovrligotey=g)is.c−caTal1hl.ecesoH,nsoiosmwtra,penmvlteeos(rtΛr,edi)naawrcskciptuihetrneaaetenorlgfyey,iqtasuctaavntlieodorwiyn- ecaRffhnRadenclgtiiemvesRi→ggrn∞av(fniRtoat=giohno0aslt(sGc)]oe;nn(scet)raanaltf,tRe(rGeliaenfftfliva=ittiyoGRniN,sl/ir(me1cR+o→v∞efRref))(dRda)o/teRsea=nroly0t
6 DE DE DE times) ; (d) |f | is small at recent times, to satisfy solar system
good accordance with current observations, Λ has some theo- R
1 and galactic scale constraints. In addition to these conditions,
retical difficulties such as its tiny value as compared with the-
: therearesomedesirablecharacteristicsthataviablecosmologi-
voretical predictions of the vacuum energy density, the cosmic
calmodelhastosatisfy(Amendolaetal.2007).Itshouldhavea
Xicoincidenceproblem,andrelatedfine-tuning.Thissituationhas radiation-dominatederaatearlytimesandasaddle-pointmatter-
motivatedthesearchforalternativeslikemodified-gravitytheo-
r dominated phase followed by an accelerated expansion as a fi-
aries.Thesimplestmodified-gravitycandidatesaretheso-called nal attractor. By using the parameters m¯ ≡ Rf /(1+ f ) and
f(R)-theories,inwhichtheLagrangiandensityL = R+ f(R)is r ≡ −R(1+ f )/(R+ f), it can be shown that,RaRn earlyRmatter-
anonlinearfunctionoftheRicciscalarR. R
dominatedepochoftheUniversecanbeachievedifm¯(r≈−1)≈
Asiswellknown,metric f(R)-theoriescanbethoughtofasa 0+ and m¯/r(r ≈ −1) > −1. Furthermore, a necessary condition
specialcaseofascalar-tensortheory;aBrans-Dickemodelwith
foralate-timeacceleratedattractoris0<m¯(r≈−2)≤1.
acouplingconstantω =0.Anacceleratedexpansionappears
BD Thereareviable f(R)gravitytheoriesthatsatisfyallthecri-
naturallyinthesetheories.Theveryfirstinflationarymodel,pro-
teriamentionedabove(Hu&Sawicki2007;Starobinsky2007;
posedbyStarobinskymorethanthreedecadesago(Starobinsky
Appleby & Battye 2007; Cognola et al. 2008; Linder 2009;
O’Dwyeretal.2013).However,thereisagenericdifficultyfrom
(cid:63) [email protected]
(cid:63)(cid:63) [email protected] which all these “viable” f(R) theories (Thongkool et al. 2009)
(cid:63)(cid:63)(cid:63) [email protected] suffer: the curvature singularity in cosmic evolution at a finite
(cid:63)(cid:63)(cid:63)(cid:63) [email protected] redshift (Frolov 2008). It can be shown that this type of sin-
Articlenumber,page1of10
A&Aproofs:manuscriptno.main
gularity problem can be cured, for instance, by adding a high- 2. γgravityreview
curvature term proportional to R2 (Appleby et al. 2010) to the
Weinvestigatespatiallyflatcosmologicalmodelsinthecontext
density Lagrangian. Therefore, it is not possible to have cos-
ofγgravity(O’Dwyeretal.2013),aviable f(R)theorydefined
micaccelerationwithatotallyconsistent f(R)theorymodifying
bythefollowingansatz:
gravity only at low curvatures. We remark that we do not ad-
dress this problem here and only consider modifications at low αR (cid:34)1 (cid:32) R (cid:33)n(cid:35)
curvatures. f(R)=− ∗γ , , (1)
n n R
Weconsiderthespecificcaseofaviable f(R)theorycalled ∗
γabglerafv(iRty)t(hOe’oDriwesy,esrtreutcatul.re2f0o1r3m).atGioenneimricpaolsley,sisnucahlmstorostngalclovni-- where γ(n,x) = (cid:82) xtn−1e−tdt is the incomplete Γ-function and
0
straintsontheparametersofthemodelsthattheeffectiveequa- α, n and R∗ are free positive constants. In reality, γ gravity can
tion of state parameter cannot be distinguished from that of a be thought of as a simple generalization of exponential gravity
cosmologicalconstant.Inγgravitythesteepdependenceonthe (Linder2009)
RicciscalarRfacilitatestheagreementwithstructureformation.
O’Dwyer et al. (2013) showed that, in principle, the parameter f(R)=−αR∗(1−e−R/R∗), (2)
thatcontrolsthesteepnessinγgravityallowsmeasurabledevia-
tionsfromΛCDM(Cosmologicalconstant,ColdDarkMatter)at obtained by fixing n = 1 in Eq.(1). We emphasize that γ grav-
ity can satisfy all the stability and viability conditions. As dis-
bothlinearperturbationandbackgroundlevels,whilestillcom-
cussedinO’Dwyeretal.(2013),forfixedn,thereisaminimum
patiblewithbothcurrentobservations.Themaingoalofthispa-
peristostudytheeffectsofγgravityonthestructureformation value(αmin)oftheparameterαsuchthatforvaluesα > αmin a
late-timeacceleratedattractorisachieved.Weconsiderthiscase
at nonlinear scales for choices of parameters where the model
hasobservablesignatures1onthebackgroundexpansionhistory throughout.FromEq.(1)weobtainthefollowingderivatives:
ofourUniverse.
Wegoonestepfurtherandanalyzethenonlinearevolutionof fR =−αe−(cid:16)RR∗(cid:17)n, (3)
sthtriuscptuapreesrcisombapsuetdedonfriosmansluigmhetrmicaoldsifiimcautliaotnioonfs.ITSIhSec(Lodlienathreast fRR = αRn(cid:32)RR (cid:33)ne−(cid:16)RR∗(cid:17)n. (4)
∗
et al. 2014), which in turn is a modification of the RAMSES hy-
drodynamic N-body code (Teyssier 2002). See also Zhao et al. WenotefromEq.(3)thatwithincreasingn,thesteepnessofthe
(2011);Lietal.(2012);Oyaizuetal.(2008);Motaetal.(2008); f(R) function increases. Higher n means smaller |f |, and the
R0
Li et al. (2011); Boehmer et al. (2010); Zumalacarregui et al. departuresfromGRwillbesmalleraccordingly.
(2013);Puchweinetal.(2013)forothercodesthathaveimple- Althoughthereisnocosmologicalconstant, f(0)=0,itfol-
mentedandperformedsimulationsof f(R)gravityand(Winther lowsfromEq.(1)thatGRwithΛisrecoveredathighcurvatures.
etal.2015)forarecentcodecomparisonbetweenthesecodes. Therefore,forR(cid:29)R themodelsbehavelikeΛCDM.Sincewe
∗
We only consider the modifications made to implement the aremainlyinterestedinphenomenathatoccurredafterthebegin-
γ gravity field equations in the modified gravity part of the N- ningofthematter-dominatedera,weneglectradiationandwrite
bodycode.Formoredetailsontheimplementationofthescalar the effective cosmological constant (the cosmological constant
fieldsandothertechnicalitieswerefertoLlinaresetal.(2014). ofthereferenceΛCDMmodel)as
For this purpose, we focus on simple observables such as the
matterpowerspectrum,halomassfunction,densityprofiles,and Λ˜ = αR∗Γ(1/n)=3H˜2(1−Ω˜ ). (5)
velocity profiles to investigate modified gravity signatures that 2n 0 m0
were previously studied in Li et al. (2012); Hammami et al.
(2015); Schmidt et al. (2009); Gronke et al. (2015); Lombriser In the equation above, Ω˜m0 denotes the present value of the
etal.(2012,2014);Wintheretal.(2012);Terukina&Yamamoto matter density parameter that a ΛCDM model would have if it
(2012); Shi et al. (2015); Pujol et al. (2014); Li et al. (2013); hadthesamematterdensitytoday(ρ¯m0)asthemodifiedgravity
Hellwingetal.(2013);Heetal.(2013,2015);Schmidt(2010); f(R)model. H˜0 representstheHubbleconstantinthereference
Braxetal.(2013);Tessoreetal.(2015);Gronkeetal.(2015);Li ΛCDMmodel.Therefore,wehavem2 ≡8πGρ¯m0/3=Ω˜m0H˜02 =
etal.(2013)thatcanalsobeobserved(Bertietal.2015;Zivick Ω H2,whereΩ andH arethepresentvalueofthematteren-
m0 0 m0 0
et al. 2015; Mak et al. 2012; Cai et al. 2014, 2015; Song et al. ergydensityparameterandHubbleparameterinthe f(R)model,
2015; Wilcox et al. 2015; Jain & Khoury 2010; Schmidt 2010; respectively.ItisusefultorewriteR as
∗
Hellwingetal.2014).
R 6nd
Thispaperisorganizedasfollows:inSect.2werevisittheγ ∗ = , (6)
gravity model and show the main properties of the background m2 αΓ(1/n)
andlinearperturbationevolution.Section3detailsthedynamics
equationsofscalaronfield(fR)andparticlemovementequations whered =(1−Ω˜m0)/Ω˜m0.Tocomputethebackgroundevolution,
forγgravity,whichmustbesolvedbyourcodeduringthesim- westartfromthe f(R)fieldequationforaFLRWmetric
ulations. The method for solving these equations is briefly ex-
Rf − f
plained in Sect. 4, and the code is tested in Sect. 5. Finally our H2(1+ f +R(cid:48)f )− R =m2e−3y, (7)
R RR
resultsareshowninSect.7,andweconcludeinSect.8. 6
where(cid:48) ≡ d/dy(y = lna), H ≡ a˙/aistheHubbleparameter(a
1 Withobservablesignatureswemeanthattheequationofstatediffers dotdenotesthederivativewithrespecttocosmictime),whichis
enoughfromtheΛCDMvalueofw = −1atlowredshiftsthatsucha relatedtoRby
deviation could be detected by near-future experiments like WFIRST
Spergeletal.(2015),forexample. R=12H2+6HH(cid:48). (8)
Articlenumber,page2of10
M.VargasdosSantosetal.:N-bodySimulationsofγGravity
Ω˜ n α
m0
0.267 2 1.05
0.267 2 1.18
0.267 2 1.5
Table1.Overviewofthemodelparametersforγgravity.
Fig.2.Fractionaldifferenceinthematterpowerspectrumwithrespect
toΛCDMfordifferentvaluesofnandα,asindicatedinTable1.These
resultsareusedtocomparewiththenonlinearpowerspectrumfromour
numericalsimulations.
In GR, f = Q = 0, and there is no scale dependence for
RR
thedensitycontrastinthelinearregime.ForΛCDMthegrowing
modecanbeexpressedintermsofhypergeometricfunction F
Fig.1.Effectiveequation-of-stateparameterw asafunctionofred- 2 1
DE as(Silveira&Waga1994)
shiftzfortheparametersgiveninTable1.Thestrongestdeviationfrom
−1islowerthan4%.
(cid:34) (cid:35)
1 11
δ+ ∝e−y2F1 ,1, ,−e3yd . (16)
Tosolvetheseequations,weintroducethenewvariables 3 6
H2 We solved Eq. (14) numerically and obtained the growing
x (y)= −e−3y−d, (9)
1 m2 mode for the γ gravity. By using (16), we then obtained the
R fractional change in the matter power spectrum P(k) relative to
x2(y)= m2 −3e−3y−12(d+x1). (10) ΛCDM. Figure 2 shows ∆Pk/PΛ at z = 0 for the three choices
ofparametersshowninTable1.
Withthesedefinitionsweobtain
x
x(cid:48)(y)= 2, (11) 3. N-bodyequations
1 3
R(cid:48) f(R)modelsareequivalenttoascalar-tensortheory(Braxetal.
x2(cid:48)(y)= m2 +9e−3y−4x2, (12) 2008), where the first derivative of the f(R) function, fR. This
fieldpropagatesaccordingtheequation
whereR(cid:48) isgivenbyEq.(7).Itisstraightforwardtoverifythat,
asdefined, x and x arealwayszeroduringtheΛCDMphase.
Weherefocu1sonthe2threecasessummarizedinTable1. (cid:3)f = ∂Veff = (1− fR)R+2f +κ2T, (17)
Furthermore, in terms of x1(y) and x2(y), the effective dark R ∂fR 3
energyequationofstate(w )isgivenby,
DE
where κ2 = 8πG/c4 and T is the trace of energy-momentum
w =−1− 1 x2 . (13) tensor, T = gµνTµν. In the quasi-static limit (see, e.g., Noller
DE 9x +d etal.(2014);Boseetal.(2015);Llinares&Mota(2014,2013))
1
thisequationbecomes
Fortheconsideredmodels,theevolutionofw asafunctionof
DE
theredshiftzisshowninFig.1. 1 R−R(a)
∇2f = −m2a−3δ , (18)
For a general f(R) model the differential equation for the a2 R 3 m
matterdensitycontrast(δ )inthelinearregimeforsubhorizon
m
scalesisgivenby(Pogosian&Silvestri2008;Zhang2006;dela where
Cruz-Dombrizetal.2008)
R(a)
=3(a−3+4d)+∆ (a), (19)
δ(cid:48)(cid:48)+(cid:32)2+ H(cid:48)(cid:33)δ(cid:48) − 1−2Q 3H02Ω˜m0 e−3yδ =0, (14) m2 R
m H m 2−3QH2(1+ f ) m and ∆ (a) = x (a)+12x (a). The Ricciscalar R in functionof
R R 2 1
f isgivenbyinvertingEq.(3)
R
where
(cid:32) α (cid:33)1/n
Q(k,y)=− 2fRRc2k2 . (15) R=R∗log |f | . (20)
(1+ f )e2y R
R
Articlenumber,page3of10
A&Aproofs:manuscriptno.main
The geodesic equation, needed to update the particle positions, 4. ImplementationintheISIScode
reads
Implementingscalar-tensortheoriesofgravityinN-bodycodeis
(cid:32) (cid:33) ratherstraightforwardbecausethescalar-tensortheoriesallcon-
1 f
x¨+2Hx˙=− ∇ Φ− R , (21) tributeasafifthforceandbecauseRAMSES,whichISISisbased
a2 2 on, has been widely used, thoroughly tested, and optimized. In
thissectionwedescribehowtheequationsweneedtosolveare
whereΦisthenewtonianpotential,whichthedynamicsisgiven implementedinISIS.Formoredetailssee(Llinaresetal.2014).
bythePoissonequation To solve for f directly is not numerically stable since the
R
solution can potentially vary over several orders of magnitude
3m2δ when going from deep voids to massive clusters in our simula-
∇2Φ= m. (22) tion. We therefore introduce a field redefinition |f | = A(u) →
2 a ∇|f | = b(u)∇u where b(u) = dA(u)/du. The Rgeneral field
R
equation for f discretized on a grid with the field-redefinition
When implementing these equations in the N-body code, we f˜ ≡ −a2f =R A(u), where u is the field we solve for, can be
R R
needtorewritethemincode-unitsgivenby
writtenas
Φa2 L(ui,j,k)=∇code·[b(u)∇codeu]i,j,k+Ωm0(H0B0)2× (31)
x˜= x/B0, Φ˜ = (H0B0)2, ×a4R∗ log(cid:34) αa2 (cid:35)1/n−a(δ ) −a−a4(cid:34)4d+ ∆R(a)(cid:35),
dt˜= H0dt, ∇ = B .∇. (23) 3m2 A(ui,j,k) m i,j,k 3 
a2 code 0
where
HereB0isthesizeofthesimulationbox.Intermsof f˜R =−a2fR, ∇ ·[b(u)∇ u] = (32)
code code i,j,k
theevolutionequationsbecomes
bi+1/2,j,k(ui+1,j,k−ui,j,k)−bi−1/2,j,k(ui,j,k−ui−1,j,k)+
h2
d2x˜ 1
dt˜2 =−∇codeΦ˜ − 2(B H )2∇codef˜R, (24) + bi,j+1/2,k(ui,j+1,k−ui,j,k)−bi,j−1/2,k(ui,j,k−ui,j−1,k)+
0 0 h2
3
∇2codeΦ˜ = 2Ωm0aδm, (25) + bi,j,k+1/2(ui,j,k+1−ui,j,k)−bi,j,k−1/2(ui,j,k−ui,j,k−1),
h2
∇2 f˜ =Ω (H B )2a4× (26)
code R m0 0 0
×−3Rm∗2 log(cid:32)αfa˜2(cid:33)1/n+(cid:34)a−3+4d+ ∆R3(a)(cid:35)+a−3δm. wanhderwehherisetwheeghraivdesdpeaficinnegdabn(du)b≡i±1/d2Ad,(uju,k).≡ 21(b(ui±1,j,k)+b(ui,j,k))
R TheequationsaresolvedusingNewton-Gauss-Seidelrelax-
ation(withmultigridacceleration).Themethodconsistsofgoing
Thesearetheonlyequationsweneedtoimplementandsolvein throughthegridandupdatingthesolutionusing
theN-bodycode.
L(u )
For comparison we also need the linearized field equation. unew =u − i,j,k . (33)
Simulationswiththisequationcomparedtothefull fR equation i,j,k i,j,k ∂L(ui,j,k)/∂ui,j,k
isagoodmeasureoftheamountofscreeningthattakesplacein
Forthiswealsoneed∂L(u )/∂u ,whichisgivenby
themodel.Thelinearized f equationissimply i,j,k i,j,k
R
∂L(u ) ∂∇ ·[b(u)∇ u]
i,j,k = code code i,j,k+ (34)
1 ∇2δf =m2(a)δf −m2a−3δ , (27) ∂ui,j,k ∂ui,j,k
a2 R φ R m −Ω (H B )2a4R∗ b(ui,j,k) log(cid:34) α (cid:35)1/n−1,
whereδf = f − f (a)andm2(a)= 1 .Incodeunits,taking m0 0 0 3m2 A(ui,j,k) A(ui,j,k) 
R R R φ 3fRR(a)
u=− δfRa2 ,weobtain where
2(H0B0)2
∂∇ ·[b(u)∇ u]
Ω a code code i,j,k = (35)
∇2 u=[m (a)aB ]2u+δ m0 , (28) ∂ui,j,k
code φ 0 m 2 − bi+1/2,j,k+bi−1/2,j,k − bi,j+1/2,k−bi,j−1/2,k+
h2 h2
andthegeodesicequationbecomes − bi,j,k+1/2+bi,j,k−1/2 + 1c (cid:34)ui+1,j,k+ui−1,j,k−2ui,j,k+
h2 2 i,j,k h2
d2x˜
dt˜2 =−∇codeΦ˜ −∇codeu. (29) +ui,j+1,k+ui,j−1,k−2ui,j,k + ui,j,k+1+ui,j,k−1−2ui,j,k(cid:35),
h2 h2
Wehave
and c = c(u ), where c(u) = db(u). Some problems arise
i,j,k i,j,k du
a2(H B )2R(a)(cid:34) R (cid:35)n related to how boundary conditions are handled on the refined
m2(a)a2B2 = 0 0 ∗ e[R(a)/R∗]n. (30) gridsinthecode.ThisisdiscussedinLlinaresetal.(2014);Li
φ 0 3αn H2 R(a)
0 etal.(2012).
Articlenumber,page4of10
M.VargasdosSantosetal.:N-bodySimulationsofγGravity
Our solver needs a starting guess, and for this we use the whereδ= ρin −1characterizesthedensitycontrastbetweenthe
cosmologicalbackgroundsolution,thatis,thesolutionwewould insideandoρouuttsideofthesphere,ρ¯ isthemeandensity,Risthe
expectwhentherearenomattersources, radiusofthesphere,and B thesizeofthebox.Thevalueofδ
0
chosen for the test is 5000. For the f(R) test we used R = 25
|fR(a)|=αe−[R(a)/R∗]n →u= A−1(|fR(a)|). (36) Mpc/hand B0 = 250Mpc/h.Asphericalsymmetricconfigura-
tion is effectively one-dimensional, therefore the field equation
This is only needed when the simulation is stared as otherwise
reducestoanODE,whichwesolvedbyusingMathematicaand
wecanusetheoldsolutionasourguess.Forγgravityourchoice
usedthistocomparewith.
forAandtherelatedexpressionsforbare
Forthesecondtestweanalyticallycomputedadensityfield
ρ(x,y,z) ≡ ρ(x) (i.e., a 1D configuration), using the field equa-
A(u)≡ f˜R(u)=αa2e−eu, (37) tion,sothatthesolutionisgivenbyasine:u∝2+sin(2πx).
b(u)= dA(u) =−αa2e−eueu, (38) Figure 3 shows the result of both these tests and the differ-
du entcolorsdepictthedifferentrefinementlevels.Weseethatthe
db(u) curvesaresmoothwhengoingfromoneleveltoanother,which
c(u)= =αa2e−eu(eu−1)eu, (39)
du demonstratesthatboundaryconditionsarehandledproperly.The
tests were performed using the serial version of the code, and
andthebackgroundvalueofuis bothtestsgivetheexpectedresults,whichdemonstratesthatthe
  codeworksproperly.
u= A−1(f˜R(a))=nloga−3+R4d/(+3H∆2R)(a)/3, (40)
∗ 0
6. Simulations
intermsofuwehave
Theγgravitysimulationswererunusing5123 darkmatterpar-
L(u)=∇[b(u)∇u]+Ωm0(H0B0)2× (41) ticleswithaboxsizeof B0 = 250Mpc/h,andwerefinedcells
(cid:34)a4R (cid:32) ∆ (a)(cid:33)(cid:35) whenever it had more than eight particles in it. The highest re-
× ∗eu/n−aδ −a−a3 4d+ R , finement level reached in our simulation was eight. The back-
3m2 m 3 groundcosmologyusedforthesimulationwascomputedusing
∂ L(u)= ∂ ∇[b(u)∇u]+Ω (H B )2(cid:32)a4R∗eu/n(cid:33). (42) Eq.(7)withh=0.71,ΩΛ =0.733andΩm0 =0.267.
∂u ∂u m0 0 0 3m2 n ThemodelparameterswerepresentedinSect.2,withsome
plotsofvariousbackgroundquantities.Wealsocomparedthese
When the fifth force is implemented, we simply replace it with simulations with simulations of the Hu-Sawicki f(R) model
aneffectiveforceFeff thatincludestheeffectsofmodifiedgrav- fromLlinaresetal.(2014).
ity wherever the code normally works with the gravitational
force,F ,
N
7. Results
Feff = FN +Fφ. (43)
7.1. Powerspectrum
The expression for the fifth force in code-units is given in
The nonlinear matter power spectrum is an important observ-
Eq.(24),andweusethesamefive-pointstencilasRAMSESuses ableandcouldbeusedtodistinguishamongdifferentmodelsof
tocomputethegravitationalforce∇Φ˜ tocomputethefifthforce
structureformation.Asweshowedabove,γ gravitycanhavea
∇f˜R. strong effect on the growth rate of the linear perturbations. We
expect these signatures to be detectable in the nonlinear matter
powerspectrum.
5. TestsoftheN-bodysolver
To compute the power spectrum we used a public code,
To verify that the solver is implemented correctly, we tested it POWMES (Colombi & Novikov2011), which uses folding meth-
severaltimes.Wepresentsomeofthesetestsinthissection.The odstocomputethepowerspectrum.
density is provided to the code through a distribution of parti- Figure 4 displays the difference of the matter power
cles. The density estimation (CIC) and refinement criteria are spectrum with respect to ΛCDM, defined as ∆P/PΛCDM ≡
thesameasthoseusedforthecosmologicalsimulations.Wealso P(k)/PΛCDM(k)−1forγgravitytogetherwiththecorresponding
have the option to set the density field analytically in the code. predictionsfromlinearperturbationstheoryandtheHu-Sawicki
To test that the treatment of the boundary of the refinement is f(R)modelforcomparison.Wefocusonthepresent-dayepoch,
correct, we included two levels of refinement. The model pa- whichcorrespondstoz=0.
rametersusedforthetestsarethesameasmentionedpreviously Figure 4 shows that the differences from ΛCDM are lower
inTable1. than5−10%forallofourrunsintherange0.05hMpc−1 (cid:46)k (cid:46)
ThefirsttestcasecorrespondstoasphereofradiusRofcon- 0.1hMpc−1)andinagreementwiththelinearperturbationresult
stantdensitylocatedinthecenteroftheboxandembeddedina seeninFig.(2).ThedeviationfromΛCDMapproacheszerofor
uniformbackground: larger scales. This is because the range of fifth force is smaller
than the horizon, therefore the modifications of gravity are not
ρ(r)= 1+(143π+δδ(cid:16))BρR0(cid:17)3, r<R , (44) feeevffleetrc,otn∆oPfth/tehPeΛlaCfirDgftMehstcfoosrnccateilneausc.etsOs,tnoansgdmrowawlelewsr,eientholankrlgiinenerγadregvsrciaaavtliieotsyn,,sti.hnHecofouwnl--l
 1+ 4πρδ(cid:16)R(cid:17)3, r>R tmernaehsctahntaocnetimhsmeenHits.u-iTnShapiwslaiiycskoyannmthioneddseiecl,aswtciaohlneerset,hwtahtheitchcheharsmecdreuelceeeonsnintgshcermepeeoncwihneagr-
3 B0
Articlenumber,page5of10
A&Aproofs:manuscriptno.main
Fig. 3. Scalar field from tests. Different colors depict the different refinement levels. The left panel shows the field from a spherical density
distribution,therightsideshowsa1Dsinefieldobtainedinthesecondtest.
Fig. 5. Fractional difference of the halo mass function for γ gravity
Fig. 4. Fractional deviation of the dark matter power spectrum with (data points) with respect to the ΛCDM model. For comparison we
respect to the ΛCDM model. The linear predictions of γ gravity are alsoshowtheresultsfromsimulationswiththeHu-Sawicky f(R)model
represented by dashed lines. For comparison we also show the re- with|f |={10−4,10−5,10−6}(dashedlines)fromLlinaresetal.(2014).
R0
sults from simulations with the Hu-Sawicky f(R) model with |f | =
R0
{10−4,10−5,10−6}fromLlinaresetal.(2014).
nismdoesnotworkverywellinoursimulation.Wediscussthis the signatures are very similar to the |f | (cid:38) 10−5 Hu-Sawicki
inmoredetailinSect.7.4. model.Forthemostmassivehalos(M/MR0 (cid:38) 1015)theeffectof
(cid:12)
screeningismuchmorepronouncedinHu-Sawickimodels,and
themass-functionapproachesΛCDM.Thisdoesnotoccurforγ
7.2. Halomassfunction
gravityandisagainanindicationthatthechameleonscreening
Thehalomassfunction,thenumberdensityofhaloswithagiven mechanismdoesnotworkveryefficientlyinoursimulations.
mass,isausefultoolforinvestigatingtheefficiencyofamodel
informinghalosofdifferentmasses.Tolocatehalosinthesimu- For γ gravity, the number of halos increases significantly,
lationoutputs,weusedtheAHF(AmigaHaloFinder(Knollmann especiallyatthehighmassend,byupto40−100%forcluster-
&Knebe2009)).Wedeterminedthehalomassfunctionbybin- sized halos. For Hu-Sawicki models, when the value of the f
R
ninghalosinlogarithmicmassintervals. field becomes comparable to the cosmological potential wells,
In Fig. 5 we show the fractional difference with respect to the chameleon effect starts to operate. This can be seen in the
ΛCDM of halo mass function computed from our simulations. massfunction,wheredeviationsfromΛCDMapproachzeroin
Our measurement of the halo mass function itself is limited by thehigh-massendformodelswith|f |=10−5and|f |=10−6.
R0 R0
statisticsandtoalesserextent,bytheresolutioninthehighand
lowmassend.However,wecanreducetheimpactofthesetwo Thelargedeviationswefind∆n/nΛCDM ∼0.5−1inthehigh-
effects by considering the relative difference between the halo massendareprobablyalreadyruledoutbypresentclustercounts
massfunctionsmeasuredinmodifiedgravityandΛCDMsimu- (Cataneoetal.2015;Maketal.2012).Thisdoesnotruleoutthe
lationswiththesameinitialcondition. model per se, but requires choices of parameters where |f | is
R0
Figure5showsthatwefindanexcessofhalosintherange muchsmallertoday,leadingtonosignaturesinthebackground
1013 (cid:46) M/M (cid:46)1015probedbyoursimulationinγgravity,and evolution(i.e.,intheHubblefactor)oftheUniverse.
(cid:12)
Articlenumber,page6of10
M.VargasdosSantosetal.:N-bodySimulationsofγGravity
Fig.6.FractionaldifferenceinthehalodensityprofileswithrespecttoΛCDMforfourdifferentmassbins.Forcomparisonwealsoshowthe
resultsfromsimulationswiththeHu-Sawicky f(R)modelwith|f |={10−4,10−5,10−6}(dashedlines)fromLlinaresetal.(2014).
R0
7.3. Haloprofiles and by 5% in the |f | = 10−6 case for all the three halo mass
R0
ranges analyzed. Only for the |f | = 10−5 parameters a mixed
R0
Wealsostudiedhowmodificationsofgravitychangethedensity
behaviorisseen,thatis,forthelowertwohalomassrangesthe
andvelocityprofilesofdarkmatterhalos.Wefocusedonhalos boostis∼ 20%,butfortheheaviesthalosthereisnodeviation
in the mass ranges 0.5−1×1013Mpc/h, 1−5×1013Mpc/h, fromΛCDMintheinnerpartsand∼ 15%highervelocitiesare
0.5−1×1014Mpc/h,and1−5×1014Mpc/h.Thedensitypro- foundintheouterparts.Forγgravity,thedifferencebetweenthe
fileswerecalculatedbybinningdarkmatterparticlesinannular
models we simulated is more expressive for less massive halos
bins for each halo. Our calculated density profiles are averages (5×1012 < M/M < 1013 and 1013 < M/M < 5×1013), for
(cid:12) (cid:12)
of all density profiles of the proper size, ranging from 10% of most massive halos (1014 < M/M < 5×1014) only the case
thevirializationradius,r =0.1Rvir,totentimesthevirialization α=1.5isdistinguishablefromthe(cid:12)others.
radius,r=10R .Thisrangewaschosentoproperlyincludeall
vir
behaviors of the fifth force on the dark matter halos while also
avoiding the inner regions of the halos, where the resolution of 7.4. Fifthforceandscreening
oursimulationsisnotsufficient.
Figure 6 shows the fractional difference with respect to Generalrelativityisverywelltestedonverysmallscales,espe-
ΛCDM in the density profiles. We first note that the inner re- cially inside the solar system. To ensure that Fφ is not relevant
atthesescales,weneedascreeningmechanismtosuppressthe
gions (R < R ) of halos for γ gravity are significantly denser
vir
than in ΛCDM. This difference is compensated for in outer re- fifthforceatsmall-scale,veryhighcurvatureregimes.
gions (R > R ). Moreover, the profiles between the different WhenγgravitywasproposedinO’Dwyeretal.(2013),the
vir
modelparametersdonotdifferappreciablyfromeachother,and authors explored the compatibility between the model and the
thispatternrepeatsforallranges.However,thedensityprofiles solar system experiments using the chameleon mechanism for
for Hu-Sawicki models in general show stronger clustering in screening,asanyother f(R)(Braxetal.2008).
thelow-densityregionsintheoutskirtsofhalosthanintheinner Asweshowedintheresultsabove,thescreeningmechanism
regions. doesnotseemtobeworkingveryefficientlyforthemodelssim-
Inthevelocitiesprofiles,showninFig.7,weexpectthefifth ulatedhere.Tocheckhowmuchthefifthforceisscreenedinour
forcetoincreasethevelocitydispersion.Thiseffectisverysim- simulations, we compared the magnitude of the fifth and New-
ilar for both models, except for most massive halos, for which tonian force on the particles in our simulation box, as in Davis
theHu-Sawickimodelsaremorescreened,causingasubstantial etal.(2012).
decreaseincomparisontoγgravity.FortheHu-Sawickimodel Figure8showsascatterplotforthiscomparisonatredshift
the velocities are boosted by ∼ 20% in the |f | = 10−4 case z = 0.Thedispersionforsmall F < 0.01isexpected(numer-
R0 N
Articlenumber,page7of10
A&Aproofs:manuscriptno.main
Fig.7.FractionaldifferenceinthevelocitydispersionprofileswithrespecttoΛCDMforfourdifferentmassbins.Forcomparisonwealsoshow
theresultsfromsimulationswiththeHu-Sawicky f(R)modelwith|f |={10−4,10−5,10−6}(dashedlines)fromLlinaresetal.(2014).
R0
ical scatter), here the forces are tiny, so the scatter here has a models considered here, Earth and Sun are almost completely
very weak effect on the simulation. The important part in this screenedandpasslocalgravityconstraintsassumingthegalaxy
figure is the behavior for large F . If we have screening, then isscreened.
N
weshouldseeF (cid:28) F /3inhigh-densityregions(whichcorre- However,fortheparametervaluesconsideredinthispaper,
φ N
spondstohighvaluesof F ).Theresultwefindshowsthatthe the screening condition gives that the galaxy is not screened
N
forceratioisroughlyonethird-whichisthelinearprediction- which again implies that the Earth and the sun is not screened
everywhere,meaningthatthereisverylittlescreeningpresentin either3.Theonlycaveattothisisthatscreeningmighthavesur-
oursimulations. vived from earlier times. At early redshift |f | (cid:28) 1 and almost
R
Tounderstandthisbetter,werevisitthescreeningconditions.
allobjectsarescreened.Theninthecosmologicalevolution|f |
Considering a spherical symmetric body with constant density R
very quickly evolves to high values |f | = O(0.01−0.1). The
ρ embeddedinthebackgroundofconstantdensityρ ,thesolu- R
c b fieldinsideourgalaxymighthavebeentrapped,ensuringscreen-
tionstothefieldequation(seee.g.Braxetal.(2012))meanthat
ing.Thisisnotexpected,buttorigorouslyruleout(orconfirm)
thefifthforceisgivenby
this possibility, we would need simulations beyond the quasi-
Fφ ∼ 13∆RRGrM2 e−mbr, sbteaytiocnadptphreoxscimopaetioofntthhiastpwaepears.sumedinoursimulations.Thisis
where
∆R |f − f |
= Rc Rb ,
R 2Φ 8. Summaryandconclusions
N
is the so-called screening factor (also called the thin-shell) and Wehaveinvestigatedthenonlinearevolutioninγgravity,a f(R)
fRb (fRc)isthevalueof fR inthebackground(insidethebody), theoryofgravitythatisaviablealternativetoΛCDM.Themod-
where ρ = ρb (ρc). If ∆RR (cid:28) 1, the fifth force is screened and els we investigated use a screening mechanism to suppress the
we recover General Relativity. If ∆R (cid:38) 1, we instead find that deviations from General Relativity at small (solar system) and
R
F ∝ 1F andgravityissignificantlymodified. large cosmological scales. Specifically, this is the chameleon
φ W3henNthefieldislocatedintheminimumofitseffectivepo- screeningmechanism.Asaresultofthisscreeningmechanism,
tential in the background2 , we can calculate the screening fac- the strongest signatures in these models are expected to occur
tor analytically. Assuming this, we find that for the γ gravity atthenonlinearregimeofstructureformation.Therefore,toun-
2 Backgroundherereferstothesurroundingsforthebodyinquestion, 3 Thisisonlytruefortheparametersconsideredinthispaper.If|f |(cid:46)
R0
notnecessarilythecosmologicalbackground. 10−5,thenthegalaxy,Earth,andourSunarescreened.
Articlenumber,page8of10
M.VargasdosSantosetal.:N-bodySimulationsofγGravity
tweenthedifferentγgravitymodelsislargerforlow-masshalos,
5×1012 < M/M <1013,andwecandistinguishthemfromeach
(cid:12)
other;thedifferencesreach∼30%closetotheboundaryregion,
R ∼ R . This difference decreases with the mass of the halos,
vir
andforthemostmassivehalos1014 < M/M <5×1014 .Only
(cid:12)
thecaseα=1.5islowerthan∼10%fromtheothers,whichare
practicallyidentical.
The chameleon mechanism - the screening mechanism that
makes f(R)gravityviable-isnotveryeffectivefortheparame-
terchoicesconsideredinthispaper.Thisexplainsthelargede-
viations from ΛCDM in the observables we considered. Espe-
ciallytheclustercountsignaturesof40−100%inthehigh-mass
endM (cid:38)1014.5M disagreewithcurrentobservations.Thisdoes
(cid:12)
notruleoutthemodelperse,butrequirechoicesofparameters
where|f |ismuchsmallertoday,whichimpliesthatthemodel
R0
has no observable signatures in the background evolution (i.e.,
Fig.8.ComparisonbetweenthemagnitudesofNewtonianforce(F ) intheHubblefactor)oftheUniverse.
N
andfifthforce(F ),thedashedlineindicatesF /F =1/3.Theforces
φ N φ
areinunitsofH /c2.
0
Acknowledgements
MVSthankstheBrazilianresearchagencyCAPESandtheUni-
veiltheimprintsofsuchtheoriesatastrophysicalscales,weran
versityofOsloforsupportandMaxB.GrönkeandAmirHam-
severalcosmologicalN-bodysimulations.Wecomparedmodels
with ΛCDM and the Hu-Sawicki model and showed that sev- mami for useful discussions. HAW is supported by BIPAC and
the Oxford Martin School. DFM would like to thank the Re-
eralastrophysicalobservables(halomassfunction,densitypro-
search Council of Norway for funding. The simulations were
files, and power spectra) show the signatures of the model sig-
performed on the NOTUR cluster HEXAGON, which is the
nificantlystrongly.
computingfacilityattheUniversityofBergen.
For the matter power spectrum we found a small deviation,
lower than 10%, on large scales (k (cid:46) 10−1 h/Mpc), which is
consistentwiththepredictionsoflinearperturbationtheory.For References
small scales (k (cid:38) 10−1 h/Mpc), on the other hand, we found a
strongincreaseinpower,thelargestdeviation(α=1.18)reaches Ade,P.A.R.etal.2014,Astron.Astrophys.,571,A22
Amendola, L., Gannouji, R., Polarski, D., & Tsujikawa, S. 2007, Phys. Rev.,
∼ 60%,andintheothercases(α = 1.05andα = 1.5)itreaches
D75,083504
∼ 50%.ThisisquitedifferentfromtheresultsfoundintheHu- Appleby,S.A.&Battye,R.A.2007,Phys.Lett.,B654,7
Sawicky model, where the screening of the fifth force is much Appleby,S.A.,Battye,R.A.,&Starobinsky,A.A.2010,JCAP,1006,005
Berti,E.,Barausse,E.,Cardoso,V.,etal.2015,ArXive-prints
more active in suppressing the enhancement of power on small
Boehmer,C.G.,Burnett,J.,Mota,D.F.,&Shaw,D.J.2010,JHEP,07,053
scales, the largest deviation is ∼ 30%, much lower than for γ Bose,S.,Hellwing,W.A.,&Li,B.2015,J.CosmologyAstropart.Phys.,2,34
gravity. Brax,P.,Davis,A.-C.,Li,B.,&Winther,H.A.2012,Phys.Rev.D,86,044015
Brax,P.,Davis,A.-C.,Li,B.,Winther,H.A.,&Zhao,G.-B.2013,J.Cosmology
Forthehalomassfunctionwefoundthatallofourrunsgave
Astropart.Phys.,4,29
verysimilarresultinthemassrange11.5<log(M/M )<14.7,
(cid:12) Brax,P.,vandeBruck,C.,Davis,A.-C.,&Shaw,D.J.2008,Phys.Rev.D,78,
and these results are very similar to what is found in the Hu- 104021
Sawickymodelfor|f |≥10−5;a10−60%increaseinthehalo Cai,Y.-C.,Padilla,N.,&Li,B.2014,ArXive-prints
R0
abundance. However, for γ gravity we found an enhancement Cai,Y.-C.,Padilla,N.,&Li,B.2015,MNRAS,451,1036
Capozziello,S.,Cardone,V.F.,Carloni,S.,&Troisi,A.2003,Int.J.Mod.Phys.,
of 40−100% in the abundance for log(M/M ) > 14.5, which
(cid:12) D12,1969
canbecomparedtotheHu-Sawickimodel,wherethehalomass Carroll,S.M.,Duvvuri,V.,Trodden,M.,&Turner,M.S.2004,Phys.Rev.,D70,
functionapproachesthatofΛCDM. 043528
For halo density profiles we found that the different runs Cataneo,M.,Rapetti,D.,Schmidt,F.,etal.2015,Phys.Rev.D,92,044009
Cognola,G.,Elizalde,E.,Nojiri,S.,etal.2008,Phys.Rev.,D77,046009
with γ gravity were not distinguished well in any mass range; Colombi,S.&Novikov,D.2011,POWMES:MeasuringthePowerSpectrumin
inallsituationswesawthattheinnerregions(R<R )ofhalos anN-bodySimulation,AstrophysicsSourceCodeLibrary
vir
for γ gravity are significantly thicker than the halos of ΛCDM, Davis,A.-C.,Li,B.,Mota,D.F.,&Winther,H.A.2012,Astrophys.J.,748,61
around10%forlessmassivehalos,5×1012 < M/M <1013,in- delaCruz-Dombriz,A.,Dobado,A.,&Maroto,A.L.2008,Phys.Rev.,D77,
(cid:12) 123515
creasing with the halo mass until they reached ∼ 30% for the Frolov,A.V.2008,Phys.Rev.Lett.,101,061103
most massive halos, 1014 < M/M < 5 × 1014. This differ- Gronke,M.,Llinares,C.,Mota,D.F.,&Winther,H.A.2015,Mon.Not.Roy.
(cid:12)
enceiscompensatedforintheouterregions(R > R ).Forthe Astron.Soc.,449,2837
vir
Hu-Sawicki model the same effect acts for less massive halos Gronke,M.,Llinares,C.,Mota,D.F.,&Winther,H.A.2015,MNRAS,449,
2837
5×1012 < M/M(cid:12) <1013,butitchangesformoremassivehalos, Hammami,A.,Llinares,C.,Mota,D.F.,&Winther,H.A.2015,Mon.Not.Roy.
5×1013 < M/M(cid:12) <1014and1014 < M/M(cid:12) <5×1014whenthe Astron.Soc.,449,3635
screeningmechanismsuppressesthefifthforce.Thisthenleads He,J.-h.,Li,B.,&Hawken,A.J.2015,ArXive-prints
He,J.-h.,Li,B.,&Jing,Y.P.2013,Phys.Rev.D,88,103507
to a suppression in the accretion of mass concentration in the
Hellwing,W.A.,Barreira,A.,Frenk,C.S.,Li,B.,&Cole,S.2014,Physical
innerregions(comparedtothenon-screenedcase). ReviewLetters,112,221102
Finally, we computed the velocity profiles, and as we ex- Hellwing,W.A.,Cautun,M.,Knebe,A.,Juszkiewicz,R.,&Knollmann,S.2013,
pected, the velocities are significantly enhanced in γ gravity JCAP,1310,012
comparedtoΛCDM.Theothercasesarecomparableinallmass Hu,W.&Sawicki,I.2007,Phys.Rev.,D76,064004
Jain,B.&Khoury,J.2010,AnnalsPhys.,325,1479
ranges, but with some important signatures. The difference be- Knollmann,S.R.&Knebe,A.2009,Astrophys.J.Suppl.,182,608
Articlenumber,page9of10
A&Aproofs:manuscriptno.main
Li,B.,Hellwing,W.A.,Koyama,K.,etal.2013,MNRAS,428,743
Li,B.,Hellwing,W.A.,Koyama,K.,etal.2013,Mon.Not.Roy.Astron.Soc.,
428,743
Li,B.,Mota,D.F.,&Barrow,J.D.2011,Astrophys.J.,728,109
Li,B.,Zhao,G.-B.,Teyssier,R.,&Koyama,K.2012,J.CosmologyAstropart.
Phys.,1,51
Li,B.,Zhao,G.-B.,Teyssier,R.,&Koyama,K.2012,JCAP,1201,051
Linder,E.V.2009,Phys.Rev.,D80,123528
Llinares,C.&Mota,D.2013,Phys.Rev.Lett.,110,161101
Llinares,C.&Mota,D.F.2014,Phys.Rev.,D89,084023
Llinares,C.,Mota,D.F.,&Winther,H.A.2014,Astron.Astrophys.,562,A78
Lombriser,L.,Koyama,K.,&Li,B.2014,JCAP,1403,021
Lombriser,L.,Schmidt,F.,Baldauf,T.,etal.2012,Phys.Rev.,D85,102001
Mak,D.S.Y.,Pierpaoli,E.,Schmidt,F.,&Macellari,N.2012,Phys.Rev.D,
85,123513
Mota,D.F.,Shaw,D.J.,&Silk,J.2008,Astrophys.J.,675,29
Noller,J.,vonBraun-Bates,F.,&Ferreira,P.G.2014,Phys.Rev.D,89,023521
O’Dwyer,M.,Joras,S.E.,&Waga,I.2013,Phys.Rev.,D88,063520
Oyaizu,H.,Lima,M.,&Hu,W.2008,Phys.Rev.,D78,123524
Perlmutter,S.etal.1999,Astrophys.J.,517,565
Pogosian, L. & Silvestri, A. 2008, Phys. Rev., D77, 023503, [Erratum: Phys.
Rev.D81,049901(2010)]
Puchwein,E.,Baldi,M.,&Springel,V.2013,MNRAS,436,348
Pujol,A.etal.2014,Mon.Not.Roy.Astron.Soc.,438,3205
Riess,A.G.etal.1998,Astron.J.,116,1009
Schmidt,F.2010,Phys.Rev.D,81,103002
Schmidt,F.2010,Phys.Rev.,D81,103002
Schmidt,F.,Lima,M.V.,Oyaizu,H.,&Hu,W.2009,Phys.Rev.,D79,083518
Shi,D.,Li,B.,Han,J.,Gao,L.,&Hellwing,W.A.2015,Mon.Not.Roy.Astron.
Soc.,452,3179
Silveira,V.&Waga,I.1994,Phys.Rev.,D50,4890
Song,Y.-S.,Taruya,A.,Linder,E.,etal.2015,Phys.Rev.D,92,043522
Spergel,D.,Gehrels,N.,Baltay,C.,etal.2015,ArXive-prints
Starobinsky,A.A.1980,Phys.Lett.,B91,99
Starobinsky,A.A.2007,JETPLett.,86,157
Terukina,A.&Yamamoto,K.2012,Phys.Rev.,D86,103503
Tessore,N.,Winther,H.A.,Metcalf,R.B.,Ferreira,P.G.,&Giocoli,C.2015,
J.CosmologyAstropart.Phys.,10,36
Teyssier,R.2002,Astron.Astrophys.,385,337
Thongkool,I.,Sami,M.,Gannouji,R.,&Jhingan,S.2009,Phys.Rev.,D80,
043523
Wilcox,H.,Bacon,D.,Nichol,R.C.,etal.2015,MNRAS,452,1171
Winther,H.A.,Mota,D.F.,&Li,B.2012,Astrophys.J.,756,166
Winther,H.A.,Schmidt,F.,Barreira,A.,etal.2015,ArXive-print
Zhang,P.2006,Phys.Rev.,D73,123504
Zhao,G.-B.,Li,B.,&Koyama,K.2011,Phys.Rev.D,83,044007
Zivick,P.,Sutter,P.M.,Wandelt,B.D.,Li,B.,&Lam,T.Y.2015,MNRAS,
451,4215
Zumalacarregui, M., Koivisto, T. S., & Mota, D. F. 2013, Phys. Rev., D87,
083010
Articlenumber,page10of10