ebook img

Suppression of large-scale perturbations by stiff solid PDF

0.18 MB·English
Save to my drive
Quick download
Download
Most books are stored in the elastic cloud where traffic is expensive. For this reason, we have a limit on daily download.

Preview Suppression of large-scale perturbations by stiff solid

Suppression of large-scale perturbations by stiff solid 5 1 0 2 Vladimír Balek∗ and Matej Škovran† n Department of Theoretical Physics, Comenius University, Bratislava, Slovakia a J 8 January 29, 2015 2 ] c q - r g Abstract [ 1 Evolution of large-scale scalar perturbations in the presence of stiff solid (solid with pressure v 2 toenergydensityratio>1/3)isstudied. Ifthesoliddominatedthedynamicsoftheuniverselong 6 enough,theperturbationscouldendupsuppressedbyasmuchasseveralordersofmagnitude. To 2 7 avoidtoo steep large-anglepower spectrum of CMB, radiation must have prevailed over the solid 0 1. long enough before recombination. 0 5 1 : 1 Introduction v i X In standard cosmology, large-scale perturbations stay unchanged throughout the Friedmann ex- r a pansion that started after inflation, except for the last period before recombination when the Newtonianpotentialwassuppressed,due to the transitionfromradiationto matter,by the factor 9/10 (see, for example, [1]). The potential is not affected even by phase transitions and annihila- tionstakingplaceinthe hotuniverse,aslongasthe matterfillingthe universecanbe regardedas idealfluid. Amongalternativescenariosconsideredintheliteraturetherearesomethatrelaxthat assumption,introducingasolidcomponentoftheuniverseformedintheearlystageofFriedmann expansion [2, 3, 4, 5, 6, 7, 8]. The solid is supposed to have negative pressure to energy density ratio w; in particular, it can consist of cosmic strings (w = 1/3) or domain walls (w = 2/3). − − Such matter starts to influence the dynamics of the universe at late times only and has no effect on the evolution of perturbations during the hot universe period. ∗e-mail address: [email protected] †e-mail address: [email protected] 1 To obtainlarge-scaleperturbations whose magnitude at recombinationdiffers from their mag- nitude at the end of inflation, we need a solid with w 1/3. A scenario with radiation-like solid ≥ (w = 1/3) was considered in [9], where it was shown that the solid produces an additional term in the gravitational potentials that can be large at the beginning but decays afterwards. If one introducesstiff solid (w >1/3)instead,the characterof the expansionofthe universechangesfor alimitedperiodandaquestionariseswhetherthis cannotcauseashiftinthe nondecayingpartof the potentials, in analogyto what we observein a universe filled with ideal fluid as it passes from one expansion regime to another due to a jump in w. If so, the incorporationof the solid into the theory,withthevalueofitsshearmodulusleftfree,wouldenlargetheintervalofadmissiblevalues ofthe primordial potential,extending insucha waythe parameterspaceofinflationaryscenarios. Apossiblerealizationofstiffsolidwouldbeasystemofequallychargedparticleswithanisotropic short-rangeinteraction. By usingYukawapotential, oneobtains stifffluid [10,11]; however,if the potential is squeezed in some direction and the particles are arranged into a lattice, the system acquires nonzero transversalas well as longitudinal sheer modulus with respect to that direction. In order that a solid, radiation-like or stiff, has an effect on large-scale perturbations, the solidification has to be anisotropic, producing a solid with flat internal geometry and nonzero shear stress. Such solidification might possibly take place in case the Friedmann expansion was precededbysolidinflation,drivenbyasolidwithw <0ratherthanbyascalarfield[12,13,14,15]. In the paper we study how a stiff solid formed during Friedmann expansion would influence theevolutionoflarge-scaleperturbations. Insection2wederivesolutionforsuchperturbationsin a one-component universe and establish matching conditions in a universe whose matter content has changed abruptly; in section 3 we determine the behavior of perturbations after the solid has been formed and find both nondecaying and decaying part of Newtonian potential after radiation prevailedagain;andinsection4wediscusstheresults. Signatureofthemetrictensoris(+ ) −−− and a system of units is used in which c=16πG=1. 2 Perturbations in the presence of solid 2.1 Evolution equations Consider a flat FRWL universe filled with an elastic medium, fluid or solid, with energydensity ρ and pressure p, and denote the conformal time by η and the scale parameter by a. Expansion of the universe is described by the equations 1 1/2 a′ = ρa4 , ρ′ = 3 ρ , (1) + (cid:16)6 (cid:17) − H where the prime denotes differentiation with respect to η, =a′/a and ρ =ρ+p. + H 2 In a perturbed universe,spacetime metric and stress-energytensoracquire smallspace-depen- dent corrections δg and δT . We will use the proper-time gauge in which δg = 0 (the µν µν 00 cosmologicaltimet= adη coincideswiththepropertimeoflocalobservers). Themetricinthis Z gauge is ds2 =ˆ a2[dη2+2B dηdxi (δ 2ψδ 2E )dxidxj], (2) ,i ij ij ,ij − − − wherethe effectiveequalityindicatesthatonlythescalarpartofthe quantityinquestionisgiven. Suppose the matter filling the universe has Euclidean internal geometry and contains no entropy perturbations. The perturbation to T is then given solely by the perturbation to g and the µν µν shift vector of matter ξ. We will use the remaining gauge freedom to impose the condition ξ =0, sothatourgaugewillbealsocomoving. Inthisgauge,theperturbationofmassdensityδρ=δT 0, 0 the energy flux density Si = T 0 and the perturbation of stress tensor δτij =δT j are [2] i i − δρ=ρ (3ψ+ ), Si =ˆ ρ B , δτij =ˆ K(3ψ+ )δ 2µET . (3) + E − + ,i − E ij − ,ij where K is the compressional modulus, µ is the shear modulus and the index ‘T’ denotes the traceless part of the matrix. (Our K is 2 times greater and our µ is 4 times greater than K and µ in [2]. We have defined them so in order to be consistent with the standard definitions in Newtonian elasticity.) The proper-time gaugeis not defined uniquely since one canshift the cosmologicaltime by an arbitrary function δt(x). Under such shift, E stays unaltered and B and ψ transform as B B+δη, ψ ψ δη, → → −H where δη =a−1δt. This suggests that we represent B and ψ as B = +χ, ψ = χ, (4) B −H where stays unaltered by the time shift and χ transforms as χ χ+δη. B → We will restrict ourselves to perturbations of the form of plane waves with the wave vector k, and eik·x. The action of the Laplacian then reduces to the multiplication by k2; in B E ∝ − particular,the definition of becomes = k2E. For simplicity,we willsuppressthe factoreik·x E E − in and , as well as in other functions describing the perturbation. They will be regarded as B E functions of η only. Evolutionofscalarperturbationsisgovernedbytwodifferentialequationsoffirstorderforthe functions and , coming from equations T µ =0 and 2G =T . The equations are [16] B E i ;µ 00 00 ′ =(3c2 +α 1) +c2 , ′ = (k2+3α 2) α , (5) B S0 − HB SkE E − H B− HE where α = ρ /(2 )2 = (3/2)ρ /ρ, c is the “fluid” sound speed (sound speed of the solid with + + S0 H suppressedcontribution of shear modulus), c2 =K/ρ , and c is the longitudinal sound speed, S0 + Sk 3 c2 =c2 +(4/3)µ/ρ . The only place where the shear modulus enters equations (5) is the term Sk S0 + c2 in the equation for . SkE B Consider a one-component universe filled with a solid that has both p and µ proportional to ρ. The quantity K is then proportionalto ρ, too, since K =ρ c2 and c2 =dp/dρ. Mechanical + S0 S0 propertiesofsuchsolidaregivencompletelybytwodimensionlessconstantsw =p/ρandµ˜ =µ/ρ. Tosimplifyformulas,wewilloftenusetheconstantβ =µ/ρ =µ˜/w ,wherew =1+w,instead + + + of µ˜. For constant w and µ˜, the quantities appearing in the equations for and are all constant, B E except for the Hubble parameter that is proportionalto η−1. Explicitly, 3 4 α= w , c2 =w, c2 =w+ β w˜, =2uη−1, 2 + S0 Sk 3 ≡ H where u=1/(1+3w). With these expressions,equations for and simplify to B E ′ =u(1+9w)η−1 +w˜ , ′ = (k2+18u2w η−2) 3uw η−1 , (6) + + B B E E − B− E and after excluding , we arrive at an equation of second order for , E B ′′+2vη−1 ′+[q2 (2v b)η−2] =0, (7) B B − − B where q = √w˜k, v = u(1 3w) and b = 24u2µ˜. The equation is solved by Bessel functions − of the argument qη, multiplied by a certain power of η. We are interested only in large-scale perturbations, that is, perturbations stretched far beyond the sound horizon. Such perturbations have qη 1, hence we can skip the term q2 in the square brackets in (7) to obtain ≪ . =η(c η−m+c η−M), (8) J Y B wheretheparametersmandM aredefinedintermsoftheparametersν =v+1/2=(3/2)u(1 w) − and n=√ν2 b as m=ν n and M =ν +n. The constants are denoted c and c to remind J Y − − us that the two terms in (8) come from the Bessel functions J and Y. The function is non-oscillating for b<ν2 and oscillating for b>ν2. Solutions of the second B kind are well defined if the solid was not present in the universe from the beginning, but was formed at a finite time. Here we will restrict ourselves to the solutions of the first kind, which means that we will consider only values of the dimensionless shear stress µ˜ (3/32)(1 w)2. ≤ − An approximate expression for is obtained by inserting the approximate expression for E B into the first equation in (6). In this way we find . =cˆ η−m+cˆ η−M, (9) J Y E where cˆ and cˆ are defined in terms of c and c as cˆ = (1/w˜)(3/2 n)c and cˆ = J Y J Y J J Y − − (1/w˜)(3/2+n)c . Y − 4 2.2 Potentials Φ and Ψ Scalar perturbations we are interested in are most easily interpreted in the Newtonian gauge, in which the metric is ds2 =ˆ a2[(1+2Φ)dη2 (1 2Ψ)dx2]. (10) − − Let us express the potentials Φ and Ψ in terms of the functions and . If we perform explicitly B E the coordinate transformation from the proper-time to Newtonian gauge, we find (see equation (7.19) in [1]) Ψ= ( E′). (11) H B− For Φ we could proceed analogically, but it is simpler to use Einstein equations. If we write the scalar part of the stress tensor as a sum of pure trace and traceless part, τij =ˆ τ(1)δ +τ(2)T, ij ,ij fromequations2G =T weobtainthatthe difference ofΦ andΨis givenby the latterquantity ij ij (see equation (7.40) in [1]), 1 ∆Φ Φ Ψ= τ(2)a2. ≡ − 2 By inserting here from the third equation (3) we find ∆Φ= µa2E. (12) − We can see that in a universe filled with an ideal fluid (µ=0) the potentials Φ and Ψ coincide. AfterinsertingintotheexpressionforΨfromthesecondequationin(5)andintotheexpression for ∆Φ from the first equation in (1), we arrive at Ψ= k−2α 2(3 + ), ∆Φ=6µ˜k−2 2 . (13) − H HB E H E For the one-component universe introduced before, expressions for Ψ and ∆Φ become Ψ= 6u2w (kη)−2(6uη−1 + ), ∆Φ=24u2µ˜(kη)−2 . (14) + − B E E With and given in (8) and (9), both Φ and Ψ are linear combinations of η−2−m and η−2−M. B E For an ideal fluid m=0 and M =2ν, so that we expect the function Φ to be linear combination of η−2 and η−2ν+, where ν+ =1+ν. This is, however,not true because the coefficient in front of η−2 turns out to be zero. Thus, if we want to establish how Φ looks like for an ideal fluid, or how Φ and Ψ look like for a solid with small µ˜, we must add the next-to-leading term to the J-partof both expressions (8) and (9). The term is suppressed by the factor (qη)2, therefore the J-part of Φ for anideal fluidis constantand the J-partof Φ andΨ fora solidwithsmall µ˜ acquiresa term proportionalto η−m. For a universe filled with an ideal fluid we have . . =η(c +c η−2ν), =cˆ +cˆ η−2ν, (15) J Y J Y B E 5 wherecˆ andcˆ aredefinedintermsofc andc ascˆ = 6uc andcˆ = 3u(w /w)c . After J Y J Y J J Y + Y − − computing the additional terms in and and inserting the resulting expressions into equations B E (14), we arrive at . Φ=C +C η−2ν+, (16) J Y where C and C are defined in terms of c and c as C = 3u2(w /ν )c and C = J Y J Y J + + J Y 12u2w νq−2c . + Y 2.3 Transitions with jump in w and µ˜ Supposethefunctionsw andµ˜ changeatthegivenmomentη (“transitiontime”)from(w ,µ˜ ) η η tr I I to (w ,µ˜ ) = (w +∆w,µ˜ +∆µ˜). (We have attached the index η to the symbols w and µ˜ in II II I I order to distinguish the functions denoted by them from the values these functions assume in a particular era.) Rewrite the first equation in (5) as 3 4 ′ =c2 (3 + )+ w 1 + β , (17) B S0 HB E (cid:16)2 η+− (cid:17)HB 3 ηE where dp dw c2 = =w +ρ η. (18) S0 dρ η dρ Becauseof the jump in w , there appearsδ-function in c2 , and to accountfor it we must assume η S0 that has a jump, too. However, on the right hand side of equation (17) we then obtain an B expression of the form “θ-function×δ-function”; and if we rewrite ′ as B d d ′ = Bρ′ = 3 ρw B, η+ B dρ − H dρ on the left hand side there appears another such expression. To give meaning to the equation we must suppose that w changes from w to w within an interval of the length ∆ρ ρ , and η I II tr ≪ send ∆ρ to zero in the end. If we retain just the leading terms in equation (17) in the interval with variable w, we obtain d dw tr η w B = + E , (19) η+ dρ −(cid:16)B 3 tr(cid:17) dρ H where we have used the fact that, as seen from the second equation in (5), the function is E continuous at η =η . The solution is tr tr + E = C . B 3 w tr η+ H Denote the jump of the function at the moment η by square brackets. To determine [ ], we s B express and in terms of w and w , compute the difference and use the I II I+ II+ II I B B B − B expression for to exclude . In this way we find I B C ∆w tr [ ]= + E . (20) I B −wII+(cid:16)B 3 tr(cid:17) H 6 Note that the same formula is obtained if we assume that the functions with jump are equal to the mean of their limits from the left and from the right at the point where the jump occurs. To justify the expression for [ ], let us compute the jump in Ψ. It holds B 3 [Ψ]= k−2 2 (3 [w ]+∆w ), −2 Htr Htr η+B Etr andifwewrite [w ]=w [ ]+∆w andinsertfor[ ],wefindthat[Ψ]vanishes. Thismust η+ II+ I B B B B be so because for Ψ we have (see equation (7.40) in [1]) 1 Ψ′′+ (2Ψ′+Φ′)+(2 ′+ 2)Ψ= δτ(1), H H H −4 where the bar indicates that the quantity δτ(1) is computed in Newtonian gauge. A jump in Ψ would produce a derivative of δ-function in the first term, but no such expression with opposite sign appears in the other terms. The jump in ′ can be found from equation (17) by computing the jump of the right hand B side, with no need for the limiting procedure we have used when determining the jump in . The B result is ∆w 5 3w 4 [ ′]=4 + − II∆w+ ∆β . (21) tr tr tr B wII+H B (cid:16) 6wII+ 3 (cid:17)E 3 Scenario with stiff solid 3.1 Expansion of the universe Suppose at some moment η the hot universe underwent a phase transition during which a part s of radiation (w = 1/3) instantaneously turned into a stiff solid (w > 1/3). In a one-component universewith givenparameterw, the density of matter falls downthe faster the greaterthe value ofw. As aresult, ifthe solidacquireda substantialpartofthe energyofradiationatthe moment it was formed, it dominated the evolutionof the universe for a limited perioduntil radiationtook over again. Let us determine the function a(η) for such universe. Denote the part of the total energy that remained stored in radiation after the moment η by s ǫ. In the period with pure radiation (η < η ) the mass density was ρ = ρ (a /a)4, so that from s s s the first equation in (1) we obtain 1 1/2 a=Cη, C = ρ a4 . (22) 6 s s (cid:16) (cid:17) In the period with a mix of radiation and solid (η >η ) the mass density is s ρ=ǫρ (a /a)4+(1 ǫ)ρ (a /a)3w+ =ρ (a /a)4[ǫ+(1 ǫ)(a /a)∆], s s s s s s s − − where ∆=3w 4. As a result, the first equation in (1) transforms into + − a′ =C[ǫ+(1 ǫ)(a /a)∆]1/2. (23) s − 7 For w>1/3 the parameter ∆ is positive, therefore the second term eventually becomes less than the first term even if ǫ 1. ≪ Supposeradiationretainedlessthanonehalfofthetotalenergyatthemomentofradiation-to- solidtransition(ǫ<1/2). The subsequentexpansionofthe universecanbe dividedinto twoeras, solid-dominatedand radiation-dominated,separatedby the time η atwhich the mass densities rad of the solid and radiation were the same. The value of η is given by rad a =a (ǫ−1 1)1/∆. (24) rad s − Suppose now that the post-transitional share of energy stored in radiation was small (ǫ 1). ≪ The universe then expands by a large factor between the times η and η , s rad . a =a ǫ−1/∆ a , rad s s ≫ andcanbedescribedinagoodapproximationasifitwasfilledfirstwithpuresolidandthenwith pure radiation. Thus, equation (23) can be replaced by a′ =. C(as/a)∆/2 for η <ηrad. (25) (cid:26)√ǫC for η >η rad The solution is . (∆/2+1)as∆/2Cη˜ ∆/12+1 for η <ηrad a= , (26) (cid:26)(cid:2)√ǫCη˜˜for η >η (cid:3) rad where η˜ and η˜˜ are shifted time variables, η˜ = η η and η˜˜ = η˜ η . From the approximate ∗ ∗∗ − − expression for a we obtain rad η˜rad = 1 ǫ−∆/∆2+1ηs, (27) ∆/2+1 and by matching the solutions at η and η we find s rad ∆/2 ∆ η = η , η = η˜ , (28) ∗ s ∗∗ rad ∆/2+1 −2 Note that equation (23) solves analytically for w = 2/3 and w = 1, when ∆ = 1 and ∆ = 2. We do not give these solutions here since will not need them in what follows. The two equations in (28) can be rewritten to formulas for the ratios of shifted and unshifted times, η˜ 1 u η˜˜ ∆ u s rad 0 = = , = +1= , η ∆/2+1 u η˜ 2 u s 0 rad where u is the value of u in the radiation-dominated era, u = 1/2. These equations stay valid 0 0 also after we replace radiationby an ideal fluid with an arbitrary pressure to energy density ratio w . To demonstrate that, let us derive them from the condition of continuity of the Hubble 0 parameter. If the universe is filled in the given period with matter with the given value of w, its scale parameter depends on a suitably shifted time η˜ as a η˜2u. Thus, its Hubble parameter is ∝ =2uη˜−1 and the requirementthat is continuous at the moment when w changes fromw to I H H w is equivalent to η˜ /η˜ =u /u . II II I II I 8 3.2 Behavior of the function B We are interested in large-scale perturbations in a universe in which the parameters w and µ˜ assume values (w ,0) before η , (w,µ˜) between η and η , and (w ,0) after η . (Most of the 0 s s rad 0 rad time we willleavew free,only atthe endwe willputw =1/3.) Denote the functions describing 0 0 theperturbationbeforeη andafterη bytheindices0and1respectively,andkeepthefunctions s rad referring to the interval between η and η without index. If only the nondecaying part of the s rad perturbation (the part with constant Φ) survives at the moment η , the functions and can s 0 0 B E be replaced by their J-parts, =c η, =cˆ = 6u c . (29) 0 J0 0 J0 0 J0 B E − Forthefunctions and wehaveexpressions(8)and(9)withηreplacedbyη˜andforthefunction B E we have the first equation (15) with c and c replaced by c and c , ν replaced by ν and 1 J Y J1 Y1 0 B η replaced by η˜˜. All we need to obtain the complete description of the perturbation is to match the expressions for , and with the help of the expressions for and at the moments η 0 1 0 s B B B E E and η . rad At the moment η , the jumps in w and µ˜ are ∆w =w w ∆w and ∆µ˜ =µ˜. By using s η η s 0 s − ≡ these values and the identity = 3 , we find 0 s 0s E − H B 1 4 [ ] =0, [ ′] = ∆w β , s s 0 B B −(cid:16)2 − 3 (cid:17)E Denote x =c . Equations for the unknowns x˜=c η˜−m and y˜=c η˜−M are 0 J0 J s Y s u 3 0 x˜+y˜= x , (1 m)x˜+(1 M)y˜= 1+8u ∆w β x , (30) 0 0 0 u − − h (cid:16)8 − (cid:17)i and their solution is u 1 u 1 0 0 x˜= (M 8uβ)x , y˜= (m 8uβ)x . (31) 0 0 u 2n − − u 2n − Atthemomentη ,thejumpsinw andµ˜ are∆w = ∆wand∆β = µ˜. Byinserting rad η η rad rad − − these values into the expressions for [ ] and [ ′] we obtain B B ∆w ∆w 5 3w 4 [ ] = + Erad , [ ′] = 4 − 0∆w+ β . rad rad rad rad rad rad B w0+(cid:16)B 3 rad(cid:17) B − w0+H B −(cid:16) 6w0+ 3 (cid:17)E H Introduce the constants X˜ =c η˜−m =p−mx˜, Y˜ =c η˜−M =p−My˜, (32) J rad Y rad where p is the ratio of final and initial moments of the period during which the solid affects the dynamicsofthe universe,p=η˜ /η˜ . Equations fortheunknownsx˜˜=c andy˜˜=c η˜˜−2ν0 are rad s J1 Y1 rad u x˜˜+y˜˜= (K X˜ +K Y˜), x˜˜+(1 2ν )y˜˜=L X˜ +L Y˜, (33) J Y 0 J Y u − 0 9 where the coefficients on the right hand side are defined as 1 ∆w K = w (m+6uw) , K =ditto with m M, J + Y w0+h − 6uw˜ i → and 8u∆w m+6uw 5 3w 4 0 L =1 m + − ∆w+ µ˜ , L =ditto with m M, J Y − − w0+ w˜ (cid:16) 6w0+ 3 (cid:17) → The solution is 1 1 x˜˜= (M X˜ +M Y˜). y˜˜= (N X˜ +N Y˜) (34) J Y J Y 2ν −2ν 0 0 with the constants M and N , α=J, Y, defined in terms of the constants L and K as α α α α u u M =L (1 2ν ) K , N =L K . α α 0 α α α α − − u − u 0 0 3.3 Behavior of potentials Knowing how the function looks like, we can establish the time dependence of the Newtonian B potential Φ and the potential describing the curvature of 3-space Ψ. Before the time η , both s potentials are the same, Φ as well as Ψ = C x . Between the times η and η , the 0 0 J0 0 s rad ∼ potentials are given by the two equations in (14) with η replaced by η˜. With and inserted B E from equations (8) and (9), both Φ and Ψ become sums of terms proportional to η˜−2−m and η˜−2−M. We have already mentioned that for µ˜ = 0 the coefficient in the first term in Φ = Ψ is zero,andone easilyverifiesthatforw >1/3andµ˜ closeto zerothe firstcoefficientinboth Φand Ψ is proportional to µ˜. (After a simple algebra we find that it is proportional to m(1 4β) 8β − − and m 8uβ for Φ and Ψ respectively, with m reducing to b/(2ν)=8uw β/(1 w) in the limit + − − β 1.) The coefficients contain the constants c and c and if we use c x˜ and c y˜ with J Y J Y ≪ ∝ ∝ x˜ and y˜ given in equation (31), we find that the second coefficient is proportional to µ˜, too. (In the expression for y˜ we encounter the factor m 8uβ again.) Both coefficients contain also the − factor x Φ , therefore for η close to η we have Φ as well as Ψ µ˜(kη˜)−2Φ . As η grows, the 0 0 s 0 ∼ ∼ first correction to the term proportional to η˜−2−m, which is of order Φ , may take over while the 0 perturbationstillremainsstretchedoverthehorizon. However,inorderthatourapproximationis valid,thistermmustbenegligibleinthefirstperiodafterthemomentη . (Notethatthisdoesnot s hold for the potential Ψ just after η : it equals Φ at η , hence it is dominated by the correction s 0 s term for a short period afterwards.) As a result, µ˜ must be not too close to zero, µ˜ (kη˜ )2. s ≫ For large enough µ˜, Φ and Ψ can become much greater in absolute value not only than Φ , 0 but also than 1. The theory then seems to collapse, but it does not because, as can be checked by direct computation, kB, ψ and remain much less than 1. (A detailed discussion for ∆w = E 0 can be found in [17].) Thus, the proper-time comoving gauge which we have implemented insteadofmorecommon,andintuitivelymoreappealing,Newtoniangauge,isnotonlyconvenient 10

See more

The list of books you might like

Most books are stored in the elastic cloud where traffic is expensive. For this reason, we have a limit on daily download.