{'fraction_non_alphanumeric': 0.03942988600290455, 'fraction_numerical': 0.01510108085182948, 'mean_word_length': 4.245735859232791, 'pattern_counts': {'":': 0, '<': 16, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Adiabatic instability in coupled dark energy-dark matter models

We consider theories in which there exists a nontrivial coupling between the dark matter sector and the sector responsible for the acceleration of the universe. Such theories can possess an adiabatic regime in which the quintessence field always sits at the minimum of its effective potential, which is set by the local dark matter density. We show that if the coupling strength is much larger than gravitational, then the adiabatic regime is always subject to an instability. The instability, which can also be thought of as a type of Jeans instability, is characterized by a negative sound speed squared of an effective coupled dark matter/dark energy fluid, and results in the exponential growth of small scale modes. We discuss the role of the instability in specific coupled CDM and Mass Varying Neutrino (MaVaN) models of dark energy, and clarify for these theories the regimes in which the instability can be evaded due to non-adiabaticity or weak coupling.


I. INTRODUCTION
In order for our cosmological models to provide an accurate fit to all current observational data, it is necessary to postulate two dramatic augmentations beyond the minimalist assumption of baryonic matter interacting gravitationally through Einstein's equations. The first assumption is that there must exist either new gravitational dynamics or a new component of the cosmic energy budgetdark matter -that allows structure to form and accounts for weak lensing and galactic rotation curves. The second assumption is that a further dynamical modification or energy componentdark energyexists, driving late-time cosmic acceleration.
In the first case, recent results [1] have added significant support to an explanation in terms of particulate cold dark matter (CDM), rather than a modification of gravity. In the case of cosmic acceleration, however, the data remains consistent with a simple cosmological constant, with modifications to gravity, or with dark energy as the correct explanation.
A logical possibility is that the two dark sectors -dark matter and dark energy -interact with each other or with the visible sector of the theory [2,3,4,5,6,7,8,9,10,11]. In fact, a number of models have been proposed that exploit this idea to address, among other things, the coincidence problem [5,6]. Further, there exist classes of modified gravity models which may either be mapped to interacting dark energy models, or closely approximated by them over a broad range of dynamical interest [11,12].
There are concerns, however, about coupling these two rather differently behaving sectors. One concern is the presence of dynamical attractors that produce a cosmic expansion history significantly different from ΛCDM ; see, for example, Refs. [4,13,14]. Another specific example is the possibility of instabilities that are not present for the uncoupled system [15,16,17].
In this paper we perform a careful analysis of the vi-ability of coupled dark energy-dark matter models. We consider implications of such coupled theories on cosmological scales, both in terms of homogeneous background expansion and the growth of linear perturbations. We focus, in particular, on a class of models in which there exists an adiabatic regime in which the dark energy field instantaneously tracks the minimum of its effective potential, as explored by Das, Corasaniti and Khoury [8].
We show that if the coupling strength is much larger than gravitational, then the adiabatic regime is always subject to an instability. The instability is characterized by a negative sound speed squared of an effective coupled dark matter/dark energy fluid, and results in the exponential growth of small scale modes. We analyze a number of different models, and show that for these models the instability strongly constrains the region in parameter space that is compatible with observations. A short version of our results was given in the recent paper [18].
We improve on previous investigations of this instability [15,16,17] in a number of ways. First, we show that the instability occurs only for coupling strengths that are strong compared to gravitational coupling, and can be evaded at weaker couplings; this point was missed all in previous work. Second, we give a simple intuitive explanation of the instability as a type of Jeans instability. Normally, for cosmological perturbations, Hubble damping converts the exponential growth of gravitationally unstable modes into power law growth. However, here the Hubble damping is ineffective, due to the fact that the effective Newton's constant for the interaction of dark matter with itself is much larger than the Newton's constant governing the background cosmology. The result is the exponential growth of perturbations. Finally, we generalize previous treatments to allow an arbitrary coupling between the dark energy sector and the visible sector, in addition to the coupling between dark energy and dark matter.
In more detail, we consider models which are characterized by a function α c (φ) governing the interaction between CDM and a quintessence type scalar field φ, and a function α b (φ) governing the interaction between visible matter (baryons) and φ. The effective Newton's constants for the interaction of CDM with itself (cc), the interaction of CDM with baryons (cb), and the interaction of baryons with themselves (bb) are (see Appendix C) (1.1c) These Newton constants are those of the Einstein frame in the short wavelength limit, and include the effect of the scalar interaction mediated by φ. Here m p is the Planck mass. The adiabatic instability occurs when G cc ≫ G, or equivalently m p α ′ c ≫ 1, when the scalar coupling between dark matter particles is strong compared to the tensor coupling 1 . In the present day Universe, this regime is excluded by observations if we assume that perturbations to the cosmological background value of φ are in the linear regime, so that the parameters (1.1) are constants (see below). First, tests of general relativity in the Solar System constrain m p α ′ b to be small compared to unity. Second, observations of tidal disruption of satellite galaxies of the Milky Way provide the constraint 2 [10] G bc √ G cc G bb − 1 0.02. (1.2) Combining these constraints excludes the regime G cc ≫ G. Nevertheless, it is still of interest to explore the adiabatic instability, for a number of reasons. First, the argument given above assumes a single type of CDM particle. In more complicated models with two different CDM species, one coupled and one uncoupled, the constraint can be evaded (see Secs. V B and V D below). Second, the Newton's constants can evolve as a function of redshift, via the dependence of φ on redshift, so present day observations (z ≈ 0) do not exclude the occurrence of the instability at high redshifts. Third, it is useful to have multiple independent constraints on coupled models. Finally, it is not necessarily true that perturbations to the cosmological background value of φ are in the linear regime. In particular, it has been claimed that observational constraints on m p α ′ b can be evaded in chameleon models 3 , due the dependence of the local effective Newton's constants (1.1) on φ and hence on the local matter density [19,20]. The structure of the paper is as follows. In the next section, we introduce the general class of models and the framework in which we work. In section III we discuss the adiabatic regime, giving both local and nonlocal conditions for its applicability. We also derive the effective equation of state for these models in the adiabatic regime; as previously noted in Ref. [8], superacceleration is a generic feature of these models. Section IV discusses the adiabatic instability that can arise in the adiabatic regime. We discuss two complementary ways of understanding the instability, a hydrodynamic viewpoint and a Jeans instability viewpoint. We also derive the range of lengthscales over which the instability operates, which was incomplete in earlier analyses [15,16]. In section V we apply our general analysis to some well-known coupled models. We perform analytic and numerical analyses of the evolution of the background cosmology and of perturbations in both coupled CDM, chameleon [19,20] and mass varying neutrino (MaVaN) [22, 23,24] models, and identify the regimes in which they are subject to the instability. We conclude in VI by summarizing our findings and discussing their implications.
In appendix A we generalize the derivation of the instability to include baryonic matter, and in appendix B we generalize the derivation of the instability to regimes where CDM cannot be treated as a pressureless fluid but instead must be treated in terms of kinetic theory. Finally, in appendix C we derive and describe the Jeans instability viewpoint on the instability.
A note on conventions: throughout this paper we use a metric signature (-,+,+,+), and we define the reduced Planck mass by m p ≡ (8πG) −1/2 .

II. MODELS WITH A COUPLING BETWEEN DARK ENERGY AND DARK MATTER
A. General class of models We begin from the following action where g ab is the Einstein frame metric, φ is a scalar field which acts as dark energy, and Ψ j are the matter fields.
The functions α j (φ) are coupling functions that determine the strength of the coupling of the jth matter sector to the scalar field. The matter sectors include cold dark matter with coupling function α c (φ), and baryons with coupling function α b (φ). This general action encapsulates many models studied in the literature. The case of equal coupling to dark and visible matter, α j (φ) = α(φ) for all j 4 , corresponds to scalar-tensor theories of gravity which have been extensively studied, recently under the name of "coupled quintessence" [4]. This class of theories includes the f (R) modified gravity theories [25]. Several authors have considered the case α b (φ) = 0, in which the dark energy couples only to dark matter. We note that this choice yields a microphysical model for Modified-Source Gravity [26] in the adiabatic regime discussed in Sec. III below (in the approximation where one considers only a dark matter source and neglects baryons).
There is good theoretical motivation for considering nontrivial and different coupling functions α j (φ), since this is a generic prediction of string theory and of higher dimensional models. In fact, typically the moduli and dilation fields of string theory must be massive today, because for massless fields it is difficult to satisfy the observational constraints [Solar System tests and fifth force experiments that require dα b /dφ 10 −2 m −1 p , equivalence principle tests that require matter coupling functions for different matter sectors aside from dark matter to be very nearly the same] in a natural way because of loop corrections. However dynamical dark energy models require massless or nearly massless fields so one is forced to confront the naturalness issue. One possible solution, suggested by Damour and Polyakov [27], is that there is an attractor mechanism under cosmological evolution that drives the theory to be very close to general relativity, α ′ j → 0 for all j. It is possible that this mechanism does not work perfectly, and that there are residual deviations in the form of nontrivial matter couplings.
In addition, while equivalence principle tests strongly constrain differences between the coupling functions for different types of visible matter, the corresponding constraints on dark matter are much weaker, as first pointed out by Damour, Gibbons and Gundlach [2]. So there is good motivation for exploring models with various couplings to dark matter.
If the jth sector consists of a fermion of mass m, the corresponding action in Eq. (2.1) can be written as with p = 2 and q = 3. Other possibilities for p and q have 4 A more precise characterization of equal coupling is α ′ j (φ) = α ′ (φ) for all j, since any constant term in α j (φ) can be absorbed into a rescaling of all dimensionful parameters in the action S j . been explored in the literature. For example, Farrar and Peebles [28] consider the case p = 0, and Bean [6] considers the case p = q. All of these choices are equivalent 5 in the nonrelativistic limit, when the fermion's rest masses dominate their gravitational interactions. The choice we make here is motivated by the fact that it satisfies the equivalence principle within the jth sector, and also the fact that it arises naturally from higher dimensional models.
The field equations resulting from the action (2.1) are and where the prime represents a derivative with respect to φ. Here we treat the matter Ψ j as a fluid with Jordanframe densityρ j and pressurep j and with a 4-velocity u j a normalized according to For much of this paper we shall work in the approximation where we consider only the gravitational effects of the dark matter and neglect the gravitational effects of the baryons. In this approximation, the sums over j in Eqs. (2.3) and (2.4) reduce to a single term involving the dark matter. Dropping the index j for simplicity, and denoting the CDM coupling function α c (φ) simply by α(φ), the field equations become and These are are the standard equations for a scalar-tensor cosmology.
For dark matter we havep = 0, and we define a rescaled density variable ρ ≡ e 3α(φ)ρ . (2.7) With this new variable, we obtain and The scalar field equation (2.9) can also be written as where the effective potential that includes the matter coupling is and the prime means derivative with respect to φ at fixed ρ. The equations of motion for the fluid are [see Ap- ∇ a (ρu a ) = 0, (2.12) and Note that the equations of motion (2.8) -(2.13) of the theory do not depend on the coupling function α b (φ) of the scalar field to visible matter, in the approximation where we neglect the baryons. Our discussion of the instability in the following sections will be valid for arbitrary α b (φ). However we note that α b (φ) will enter when we want to compare with observations, since the metric that is measured is e 2α b g ab . When discussing observations below we will focus on the case α b = 0, where the dark energy is coupled only to dark matter.

III. THE ADIABATIC REGIME
The effective potential (2.11) governing the evolution of the scalar field is the sum of two terms, one arising from the original potential V (φ), and the other arising from the coupling to the energy density of the dark matter fluid. It is possible, for appropriate choices of the potential and the coupling function, that competition between these terms leads to a minimum of the effective potential. Further, in some regimes, it may be possible for the solution of the equation of motion for φ to adiabatically 6 in the track the position of this minimum. That is to say, the timescale or lengthscale for φ to adjust itself to the changing position of the minimum of the effective potential may be short compared to the timescale or lengthscale over which the background density is changing.
This adiabatic regime has been previously discussed for spatial variations of φ in the interior of massive bodies in the so-called Chameleon field models [19,20], and for the time variation of φ in a cosmological context by Refs. [8,16]. It also has been studied for the specific case of f (R) modified gravity models [29,30,31], for which the action in the Einstein frame is of the general form (2.1). In this section we review the adiabatic regime in a general context and give a careful derivation of its domain of validity.
The adiabatic approximation consists of (i) omitting the d'Alembertian term in Eq. (2.10), which gives an algebraic equation for φ that one can solve to obtain φ as a function of the density ρ; (ii) omitting the terms involving the gradient of φ from the field equation (2.8). The resulting equations are the same as those of Modified-Source gravity [26], in the approximation where one considers only a dark matter source and neglects baryons. The equations can be written as general relativity coupled to a fluid with an effective energy density ρ eff and effective pressure p eff : Here φ m (ρ) is the solution of the algebraic equation for φ. Eliminating ρ between Eqs. (3.2) and (3.3) gives the equation of state p eff = p eff (ρ eff ).
In the adiabatic regime, the matter and scalar field are tightly coupled together and evolve as one effective fluid. By taking the divergence of the field equation (3.1), we see that this fluid obeys the usual fluid equations of motion with the given effective equation of state. In a cosmological context, the effective fluid description (3.1) is valid for the background cosmology and for linear (and, indeed, nonlinear) perturbations. Therefore, the equation of state of perturbations is the same as that of the background cosmology, or the so-called entropy perturbation vanishes.

A. Condition for global validity of adiabatic approximation
Of course, the adiabatic approximation will not be a good one for all choices of V (φ) and α(φ), or in all physical situations. Our goal in this subsection is to establish criteria under which we can trust the adiabatic approximation. Roughly speaking, the mass of the scalar field associated with the minimum of the effective potential must be sufficiently large. More precisely, we define the effective mass as a function of density ρ by where the derivatives are taken at constant ρ, and φ m (ρ) is the value of the scalar field which minimizes the effective potential at a given density ρ. We assume that m 2 eff is positive, otherwise there is no local minimum of the potential and an adiabatic regime does not arise. For a perturbation with lengthscale or timescale L and density ρ to be in the adiabatic regime, it is necessary that We call this condition the local adiabatic condition. 7 If the condition is satisfied everywhere in spacetime, then the adiabatic approximation will be good everywhere. However if the condition is satisfied in a local region, it does not necessarily follow that the adiabatic approximation is valid in that region, due to non-local effects [19,20,31]. This is discussed further in Sec. III B below. We now derive a slightly more precise version of the local adiabatic criterion (3.6). For a given density ρ 0 we define φ 0 = φ m (ρ 0 ) and δφ = φ − φ 0 . We expand the potential as and defining W (φ) = e α(φ) we similarly expand The effective mass (3.5) is then given by 9) and the condition (3.4) that the effective potential be minimized yields Also a short computation using the definition (3.4) of the function φ m (ρ) gives where α 0 = α(φ 0 ) and α ′ 0 = α ′ (φ 0 ). We also define L to be the smallest lengthscale or timescale over which the density ρ changes, so that ∇ a ∇ a ρ ∼ ρ/L 2 and (∇ρ The field equations are the Einstein equation (2.8), the scalar field equation (2.9), and the fluid equations (2.12) and (2.13). These equations are not all independent, since the Einstein equation enforces conservation of the total stress energy tensor. We will take as the independent equations just the Einstein equation (2.8) and 7 It has also been called the Compton condition [31], since the RHS of Eq. (3.6) is the effective Compton wavelength of the field. the fluid equations (2.12) and (2.13), since the scalar field equation can be derived from these. 8 Therefore, to justify the adiabatic approximation, it is sufficient to justify dropping the scalar field derivative terms in the Einstein equation (2.8). The ratio of the scalar field gradient terms to the potential term evaluated at the adiabatic solution is of order By combining Eqs. (3.10) and (3.11) we obtain 13) and using this to eliminate one of the factors of dφ m /d ln ρ from Eq. (3.12) gives Thus, the adiabatic approximation will be valid whenever Now, since the first factor on the left hand side involves derivatives of logarithmic factors, we expect this prefactor to generically be of order unity. When this is true the condition (3.15) reduces to the local adiabatic condition (3.6).

B. Nonlocal condition for breakdown of adiabatic approximation
As mentioned above, the condition L ≫ m −1 eff in a local region is a necessary but not a sufficient condition for the validity of the adiabatic approximation in that region. This is because the corrections to the adiabatic approximation are determined by a wave equation obtained by perturbing Eq. (2.9) whose solutions depend in a non-local way on its sources. The corrections in a given region can become large due to a breakdown of the local adiabatic condition (3.6) that occurs elsewhere [31].
For the special case of static, spherically symmetric systems, and for a constant density object in a constant density background, the chameleon field papers [19,20] derived a precise condition for the validity of the adiabatic approximation inside the system, the so-called "thin shell" condition, which depends on the asymptotic value of the potential. The thin-shell condition is both a necessary and sufficient condition, but it restricted to static situations. Below we will show that for several specific examples the thin-shell condition and the adiabatic condition (3.6) give the same predictions in order of magnitude for the boundary of the adiabatic regime.
We now discuss an order-of-magnitude non-local criterion for static, spherically symmetric situations that predicts when the adiabatic approximation breaks down even when the local adiabatic condition (3.6) is satisfied. The condition is a generalization of a condition derived by Hu [31] in the context of f (R) modified gravity models, and is also a generalization of the thin-shell condition of Refs. [19,20].
Consider a spherically symmetric density profile ρ(r), which we will assume for simplicity is monotonically decreasing as r increases. From this density profile we can compute the corresponding adiabatic scalar field profile (3.16) If the local adiabatic condition (3.6) is satisfied for all r, then this adiabatic field is a good approximation to the actual solution φ(r). Suppose therefore that the local adiabatic condition is violated 9 for some set of values r 1 < r <r 2 . Then typically what occurs is that there is a larger interval r 1 < r < r 2 with r 1 <r 1 and r 2 >r 2 on which the solution φ(r) differs significantly from the adiabatic field profile (3.16). For a given interval (r 1 , r 2 ), significant deviations from the adiabatic approximation should occur in the vicinity of r = r 1 if where The criterion assumes that r 2 ≫ r 1 and ρ 2 ≪ ρ 1 .
The criterion (3.17) is a generalization of other criteria that have appeared in the literature. First, if one assumes that the density varies on a lengthscale of order ∼ r, then using this to approximate the integral in (3.17) gives is the condition derived in Ref. [31] in the context of f (R) models, although there the factor of α ′ 2 e α2 was neglected. The condition (3.18) was found to reliably predict the onset of deviations from the adiabatic profile in numerical solutions for the Solar System and the Galaxy [31].
Second, one can derive from Eq. (3.17) the thin-shell condition of Refs. [19,20], up to a factor of order unity. Suppose one has a uniform density sphere of density ρ c and radius R c , embedded in a uniform density, infinite medium of density ρ ∞ . For this case the local adiabatic criterion is satisfied everywhere except at the point of discontinuity of the density at r = R c , so we expect a breakdown of the adiabatic approximation near this point. Following Ref. [19] we assume that α ′ = β/m p is a constant, and that we are in the regime where α(φ) ≪ 1, so e α ≈ 1. If we apply the criterion (3.17) to an interval (r 1 , r 2 ) with r 1 < R c < r 2 , we obtain that the adiabatic approximation should fail near r = r 1 if where Φ c ∼ ρ c R 2 c /m 2 p is the Newtonian potential at the center of the sphere. If the left hand side is large compared to unity, then Eq. (3.19) will be satisfied for all values of r 1 except very close to r 1 = 0, and it follows that the adiabatic approximation will not apply in most of the interior of the sphere. On the other hand, if the left hand side of Eq. (3.19) is small compared to unity, then it follows that one would expect the adiabatic approximation to be valid throughout the interior of the sphere except very close to the surface, in a thin shell of thickness ∼ R c ∆φ/(βm p Φ c ). Both of these conclusions agree in order of magnitude with those of Refs. [19,20] 10 .
We now turn to the derivation of the non-local condition (3.17). We start by noting that static solutions of the equation of motion (2.10) can be obtained by extremizing the energy functional The basic idea [31] is that the adiabatic field profile (3.16) minimizes just the potential energy, and there some kinetic energy cost for following this adiabatic profile. When this kinetic energy cost becomes sufficiently large, it becomes energetically favorable for the field to switch to a different, non-adiabatic profile with a smaller kinetic energy and with a larger potential energy, for a net gain in energy. We now compute the total energies E ad for the adiabatic profile and E trial for an alternative trial profile which is qualitatively similar to numerical solutions [31]. When there is a net gain in energy for switching to the trial profile, and so the adiabatic profile is no longer a good approximation. The trial profile is simply φ trial (r) = φ 2 = constant for r 1 + ∆r ≤ r ≤ r 2 , and φ trial (r) = α − β/r for r 1 ≤ r ≤ r+∆r, with α and β chosen to satisfy continuity at r 1 and r 1 +∆r, φ trial (r 1 ) = φ 1 and φ trial (r 1 +∆r) = φ 2 . This trial profile depends on one parameter, namely the width ∆r of the region in which the field transitions from φ 1 to φ 2 . We obtain In deriving this formula we have used an approximate Taylor expansion of the effective potential about the point φ 2 together with Eq. (3.4): (3.23) Consider first the kinetic energy contributions, the first and third lines of Eq. (3.22). The minimum kinetic energy for any field configuration which satisfies φ(r 1 ) = φ 1 and φ(r 2 ) = φ 2 is the energy of the Laplace equation solution with these boundary conditions 11 , namely Therefore we can write the kinetic energy of the adiabatic field profile [the negative of the first line of Eq.
where ε is dimensionless and nonnegative. We will assume that ε 1, since ε ≪ 1 would require the adiabatic profile φ ad (r) to be very close to the 1/r profile, which would be a fine tuning. For generic density profiles we expect ε 1. If we now choose ∆r = 2r 1 /ε, then it follows that the net gain in kinetic energy is using r 2 ≫ r 1 . Turn now to the potential energy contributions, the second line of Eq. (3.22). If we assume that the potential V (φ) is monotonically decreasing, as in the Chameleon field models, then it follows that φ ad (r) is an increasing function of r, so that This implies that the potential energy term in Eq. (3.22) is bounded above by Inserting this together with the estimate (3.24) of the kinetic energy gain into Eq. (3.21) yields the criterion (3.17). 11 This can be shown using the Schwarz inequality ( R f gdr) 2 ≤ R f 2 dr R g 2 dr, with f = rφ ′ (r) and g = 1/r.

C. Observed equation of state parameter
In this section we derive the observed equation of state parameter w obs for these models, in the adiabatic regime, for the case α b = 0. Generically we have w obs < −1 corresponding to superacceleration, as previously noted in Ref. [8].
The background cosmological evolution is given in the Einstein frame by the equation where ρ 0 is a constant, together with the evolution equation for the scalar field. Observations of the acceleration of the Universe are fit to the model where ρ 1 = e α0 ρ 0 is the observed matter density today 12 , α 0 is the value of α today, and ρ DE is the inferred "dark energy density". The equation of state parameter is then given by Combining Eqs. (3.26) -(3.28) and using Eq. (3.4) to eliminate the term proportional to dφ/da gives .

(3.30)
This formula has the property that w obs = −1 today for all models. Also expanding to first order about a = 1 we obtain 1/w obs = −1 + ln

IV. ADIABATIC INSTABILITY
In the adiabatic regime, the models (2.1) discussed here can exhibit instabilities on small scales characterized by a negative sound speed squared of the effective coupled fluid. This instability extends down to the smallest scales for which the adiabatic approximation is valid. 12 Note that in this model some fraction of the mass density in the V term in Eq. (3.26) will cluster, so measurements of mass density using clusters will not correspond exactly to measurements of the second term in Eq. (3.26). We will not address this issue here.
Starting with a uniform fluid, the instability will give rise to exponential growth of small perturbations. The final state of the coupled fluid is beyond the scope of this paper. Theories that exhibit this instability are typically ruled out as models of dark energy.

A. Hydrodynamic viewpoint
This adiabatic instability was first discovered by Afshordi, Zaldarriaga and Kohri [15] in a context slightly different to that considered here. That context was the mass varying neutrino model of dark energy, where a dynamical dark energy model is obtained by coupling a light scalar field to neutrinos but not to dark matter. The instability was previously discussed in the context of the models considered here by Kaplinghat and Rajaraman [16]. In this paper we will generalize the treatment of the instability given in Ref. [16].
A fairly simple form of the stability criterion can be obtained by writing the potential V (φ) as a function V (α) of the coupling function α(φ) by eliminating φ. This gives where we have used Eqs. (3.2) and (3.4). The square c 2 a of the adiabatic sound speed is then given by (4.2) The system will be unstable if Equation (4.2) gives a simple prescription for computing when a given theory will be unstable, assuming it is in the adiabatic regime. This equation furnishes c 2 a as a function of α, which can be re-expressed as a function of φ using α = α(φ), and then as a function of density ρ using φ = φ m (ρ) from Eq. (3.4).
More generally, outside of the adiabatic regime, the effective sound speed is the ratio of the local pressure and density perturbations and can differ from c 2 a . In a cosmological context we have c 2 a ≡Ṗ /ρ and c 2 s (k, a) ≡ δP (k, a)/δρ(k, a), where k is spatial wavenumber and a is scale factor. For most of this paper we will consider only the adiabatic regime in which c 2 s → c 2 a , although the more general regime will be probed in our numerical integrations in Sec. V.
We now argue that the adiabatic sound speed (4.2) is always negative in the adiabatic regime. In the definition (3.5) of the effective mass, we rewrite the φ derivative in terms of α derivatives using α = α(φ). Using the fact that the first derivative of V eff with respect to φ vanishes at φ = φ m (ρ), we obtain where the derivatives are taken at constant ρ. Using Eq.
(2.11) and then eliminating ρ using Eq. (3.4) gives Simplifying using the expression (4.2) for the sound speed squared gives The adiabatic instability arises in the case c 2 s < 0 and m 2 eff > 0. [When m 2 eff < 0 there is no adiabatic regime, in the sense that adiabatic solutions are subject to a tachyonic instability. This instability has been discussed in the context of f (R) gravity models in Refs. [32,33,34].] It follows from Eqs. (4.2) and (4.6) that the adiabatic instability occurs when or when However, the region in the (V ,α , V ,αα ) parameter space defined by the union of the regions (4.7) and (4.8) is just the adiabatic regime; the adiabatic approximation requires both that V ,α < 0 so that the effective potential has a local extremum, and also that m 2 eff > 0 so that the local extremum is a local minimum. It follows that the adiabatic instability condition (4.3) is satisfied throughout the adiabatic regime. This genericity of the instability has also been deduced by Kaplinghat and Rajaraman [16] using a different method. However, as we discuss in the following subsection, the instability actually occurs only if the coupling α ′ (φ) is sufficiently large, a point missed in Refs. [15,16].

B. Scales over which the instability operates
For the instability to be relevant on some spatial scale L, then we must have in order that the adiabatic approximation be valid, as discussed above. There is also an upper bound on the range of spatial scales which comes about as follows. When the instability is present, spatial Fourier modes with wavelength L grow exponentially on a timescale where c 2 s is given by Eq. (4.2). If this timescale is longer than the Hubble time H −1 , then the mode does not have time to grow and the instability is not relevant. Therefore the range of scales over which the instability operates is More generally, for a fluid of density ρ, if the instability is to be unmodified by the gravitational dynamics of the fluid, then the instability timescale must be shorter than the gravitational dynamical time, which from Eq.
Here ρ eff (ρ) is the total mass density (3.2) of the coupled dark matter-dark energy fluid. This gives the criterion Here the quantity d ln V /dα(α) on the right hand side is expressed as a function of φ using α = α(φ), and then as a function of the density using φ = φ m (ρ). At longer lengthscales, it is possible that the negative sound speed squared still engenders an instability, but determining this requires a stability analysis of the system in question including the effects of self gravity. The results might vary from one system to another. In this paper we will restrict attention to the regime (4.12) where the presence of the instability can be easily diagnosed.
The factor in square brackets in Eq. (4.13) is always smaller than unity, since V is assumed to be a decreasing function of φ and hence also of α. It follows that the ratio of the maximum lengthscale L max to the minimum lengthscale L min satisfies Hence in order for there to be a nonempty regime in which the instability operates, the strength of the coupling of the field to the dark matter must be much stronger than gravitational strength, as discussed in the introduction. In Sec. IV C below we give a simple explanation for this requirement. In Sec. V D below we give an example of a model where we confirm numerically instability is present at strong coupling but not when the coupling is weak.

C. Jeans instability viewpoint
There are two different ways of describing and understanding the instability, depending on whether one thinks of the scalar-field mediated forces as being "gravitational" forces or "pressure" forces. From one point of view, that of the Einstein frame description, the instability is independent of gravity. This can be seen from the equation of motion (3.1); the instability is present even when the (Einstein-frame) metric perturbation due to the fluid can be neglected. In the adiabatic regime the acceleration due the scalar field is a gradient of a local function of the density [cf. Eq. (A10) below], which can be thought of as a pressure. The net effect of the scalar interaction is to give a contribution to the specific enthalpy h(ρ) = dp/ρ of any fluid which is independent of the composition of the fluid. If net sound speed squared of the fluid is negative, then there exists an instability in accord with our usual hydrodynamic intuition.
From another point of view, however, that of the Jordan frame description, the instability involves gravity. The gravitational force in this frame is mediated partly by a tensor interaction and partly by a scalar interaction. The effective Newton's constant describing the interaction of dark matter with itself is where k is a spatial wavevector. This is Eq. (C1) of Appendix C specialized to i = j = c, with α c written just as α, and specialized to α b = 0 (since experiments tell us that |α ′ b | 10 −2 m −1 p today). Here the 1 in the square brackets describes the tensor interaction and the second term the scalar interaction. At long lengthscales, k ≪ m eff /(m p |α ′ |), the scalar interaction is suppressed and we have G cc ≈ G. At short lengthscales, k ≫ m eff , the scalar field is effectively massless and G cc asymptotes to a constant, G cc ≈ G[1 + 2m 2 p (α ′ ) 2 ]. However, when m p |α ′ | ≫ 1 there is an intermediate range of lengthscales, in which the effective Newton's constant increases linearly with k 2 , A gravitational interaction with G cc ∝ k 2 behaves just like a (negative) pressure in the hydrodynamic equations. This explains why the the effect of the scalar interaction can be thought of as either pressure or gravity in the range of scales (4.17). Note that the range of scales (4.17) coincides with with the range (4.13) derived above, up to a logarithmic correction factor. From this second, Jordan-frame point of view, the instability is simply a Jeans instability. In a cosmological background with Einstein-frame metric ds 2 = −dt 2 +a(t) 2 dx 2 , the evolution equation for the CDM fractional density perturbation δ with comoving wavenumber k c on subhorizon scales in the adiabatic limit is [37,38] δ + 2Hδ − 4πG cc e α ρδ = 0, (4.19) where H =ȧ/a. Here G cc is given by the expression (4.16) evaluated at the physical wavenumber k = k c /a, and we have neglected photons and baryons. Now in the absence of the Hubble damping term in Eq. (4.19), the gravitational interaction described by the last term would cause an exponential growth of the mode, the usual Jeans instability of uniform fluid. Normally in a cosmological context, the Hubble damping term is present and the timescale ∼ 1/H associated with this term is of the same order as the timescale 1/ √ Gρ associated with the gravitational interaction in the last term. Because of this equality of timescales, the exponential growth is converted to power law growth by the Hubble damping. In the present context, however, things work differently. The gravitational constant governing the gravitational self-interaction of the mode is G cc (k) instead of G, and consequently the timescale associated with the last term in Eq. (4.19) is shorter than the Hubble damping time by a factor of Therefore the Hubble damping is ineffective and the Jeans instability causes approximate exponential growth rather than power law growth.
The above discussion can also be cast in terms a scaledependent sound speed instead of a scale dependent Newton's constant; this clarifies the relation to our previous discussion of Secs. IV A and IV B. As a slight generalization, consider a fluid with an intrinsic sound speed c s,in . Then the evolution equation (4.19) generalizes to (see Appendix A) where the effective total sound speed squared is In the range of lengthscales (4.17) this squared sound speed is a constant, independent of k, as for a normal, hydrodynamic sound speed. Outside of this range of scales, we have c 2 tot ∝ 1/k at both large and small k (if the intrinsic sound speed can be neglected), describing a conventional gravitational interaction.
We reiterate that the existence of the range of scales (4.17) in which Newton's constant scales linearly with k 2 depends on the assumption of strong coupling, |α ′ |m p ≪ 1. If, instead, |α ′ |m p 1, the dependence of G on k is very close to that of standard gravity, and the instability reduces to the normal Jeans instability of a fluid, characterized in a cosmological context by power law growth. Up till now we have described cold dark matter as a pressureless fluid. However, at a more fundamental level, one should use a kinetic theory description based on the collisionless Boltzmann equation. In the conventional ΛCDM framework, the fluid approximation breaks down at small scales, below the free-streaming lengthscale, and also in the nonlinear regime after violent relaxation has taken place in CDM halos [35]. We now discuss how, in the models discussed here, the conventional picture for the fluid domain of validity is slightly modified.
Let us denote by σ the rms velocity of the dark matter particles. Consider a perturbation characterized by a wavelength λ and wavenumber k = 2π/λ. The characteristic growth or oscillation time associated with this perturbation is τ (k) ∼ λ/ |c tot (k)| 2 , where the total effective sound speed c tot is given by Eq. (4.22). The distance traveled by a dark matter particle in this time is d(k) ∼ στ (k), and the ratio of this distance to the size of the perturbation is When this dimensionless ratio is of order unity or larger, perturbations do not have time to grow before they are washed out by free streaming of the particles, and the fluid approximation breaks down. Using the formula (4.22) with c s,in set to zero 13 , , we obtain (4.24) Now in the conventional CDM framework, the factor in the square brackets is unity, so the ratio (4.24) is proportional to k and becomes large as k → ∞. Hence the fluid approximation breaks down on small scales, below the critical free-streaming lengthscale λ FS ∼ σ/ √ Gρ. In the present context things work a little differently due to the scale dependence of Newton's constant. We can rewrite Eq. (4.24) in the approximate form range of scales (4.17), or equivalently the sound speed discussed in Sec. IV A. We see that the ratio d(k)/λ is proportional to k at large scales and at small scales, but that in the intermediate range of scales it is a constant, so that the effect of free streaming is equally important for all the modes in this range. If we define the free streaming lengthscale λ FS (σ) to be the smallest lengthscale for which free streaming is unimportant, d(k)/λ 1, then we obtain This lengthscale jumps discontinuously at σ 2 ∼ |c 2 s |. There are thus two different regimes that occur: • When σ 2 ≪ |c 2 s |, free streaming is important only at scales small compared to 1/m eff (for which the adiabatic approximation is invalid anyway). The fluid approximation is valid throughout the range of lengthscales (4.17), and so the adiabatic instability is present. This conclusion is confirmed by a kinetic theory analysis (see Appendix B), which shows that linearized perturbations of any homogeneous, isotropic initial particle distribution function are always unstable on scales that are in the regime (4.17), as long as σ 2 ≪ |c 2 s |.
• When σ 2 |c 2 s |, free streaming becomes important and the fluid approximation breaks down throughout the range of scales (4.17). One expects the free streaming (also called Landau damping) to kill the instability. This is confirmed by our kinetic theory analysis of Appendix B: we show that for a Maxwellian distribution, the finite velocity dispersion stabilizes the coupled fluid whenever σ 2 ≥ |c 2 s |, in agreement with the analysis of Ref. [15].
Consider now the evolution of cosmological perturbation modes. When will our analysis of instability apply? First, the condition σ 2 ≤ |c 2 s | is not very restrictive, since the CDM cools rapidly with the Universe's expansion, with temperature scaling as (1 + z) 2 once the particles become non-relativistic. For example, for CDM particles of mass ∼ GeV with weak scale cross sections, the CDM temperature is ∼ 10 −4 K at decoupling [35]. However, after perturbations go nonlinear and violent relaxation takes place in CDM halos, the effective coarse-grained velocity dispersion becomes much larger. Hence, our analysis does not apply to modes that are in the nonlinear regime. Our analysis will apply in the early Universe, before any modes have gone nonlinear. It will also apply to large scale modes, even after smaller scale modes have gone nonlinear, since such large scale modes should still be well described by linear theory (see for example the qualitative arguments in chapter 28 of Peebles [36]). Our investigation of specific models later in the paper will focus on these large scale, linear regime modes.
We note that earlier investigations of the instability focused instead on small scale modes, below the free streaming scale λ FS [15,16]. From the formula (4.26) we see that, for the models discussed here, either the adiabatic approximation is not valid on these small scales since λ FS 1/m eff , or σ 2 |c 2 s | and the instability is killed by free streaming.

V. EXAMPLES OF THEORIES WITH ADIABATIC INSTABILITY
In this section we discuss some specific classes of theories.

A. Exponential potential and constant coupling
We first consider theories with exponential potentials of the form with λ > 0 and with linear coupling functions where β = 2/3 and C is a constant. 14 These theories have been previously studied in Ref. [4]. The effective potential is, from Eq. (2.11), and solving for the local minimum of this potential yields the relation between φ and ρ in the adiabatic regime: Note that C must be negative in order for the effective potential to have a local minimum and for an adiabatic regime to exist. We will restrict attention to this case, and we define the dimensionless positive parameter γ = −λ/βC. The corresponding effective mass parameter is Next, we compute the sound speed squared. Using Eq. We now insert the effective mass (5.5), the sound speed squared (5.6) and the logarithmic derivative (5.7) into Eq. (3.15) and into the second half of Eq. (4.12). This yields the range of spatial scales L over which the instability operates for a given density ρ to be L min (ρ) ≪ L ≪ L max (ρ), where and L max (ρ) 2 is a constant times this: Thus, there is a nonempty unstable regime only when β|C| ≫ 1, ie with the scalar coupling is strong compared to the gravitational coupling, in agreement with the discussion in Sec. IV B.
To see the effect of the instability more explicitly, we consider cosmological perturbations. The Einstein-frame FRW equation in the adiabatic limit is where ρ ∝ 1/a 3 . This yields a(t) ∝ t 2/(3+3w eff ) , where the effective equation of state parameter is In the strong coupling limit |C| → ∞ that we specialize to here, w eff → −1. Thus the adiabatic regime of this model with large |C| is incompatible with observations in the matter dominated era, where we know w eff ≈ 0 except for at small redshifts. Nevertheless, the model is still useful as an illustration of the instability. We find from Eqs. (5.8), (5.9) and (5.10) that the range of unstable scales is given by where k is comoving wavenumber. This range of scales always lies just inside the horizon. A given mode k will evolve through this unstable region before it exits the horizon.
Next we use the approximate form (4.18) of Newton's constant in the perturbation evolution equation (4.19), and transform from t derivatives to a derivatives. This gives Specializing this equation to the exponential model using Eqs. (5.2), (5.4), (5.5) and (5.10) and taking the strong coupling limit |C| → ∞ gives (5.14) In the strong coupling limit H is approximately a constant, H ≈ H 0 , and the growing mode solution is where K 1 is the modified Bessel function. The mode grows by a factor ∼ e when the scale factor changes from a to a + ∆a, where ∆a/a ∼ aH 0 /k ≪ 1 for subhorizon modes.
A more detailed analysis of the cosmology of this model is given in Ref. [38], but in the non-adiabatic regime |C| ∼ 1 rather than the strong coupling regime |C| ≫ 1 considered here.

B. Two component dark matter models
We next consider models in which there are two dark matter sectors, a density ρ c which is not coupled to the scalar field, and a density ρ co which is coupled with coupling function (5.2) and exponential potential (5.1). Both of these components are treated as pressureless fluids. The FRW equation for this model in the adiabatic limit is [cf. Eq. (5.10) above] Similar two component models have been considered by Farrar and Peebles [28]. This model is also similar to the mass varying neutrino model model [15,17] where the neutrinos play the role of the coupled component; see Sec. V D below. The first two terms on the right hand side of Eq. (5.16) act like a fluid with equation of state parameter given by (5.11), and in the strong coupling limit |C| ≫ 1 this fluid acts like a cosmological constant. Thus, the background cosmology can be made close to ΛCDM by taking |C| to be large. The fraction of dark matter which is coupled must be small in the limit of large coupling, |C| ≫ 1. Denoting , Ω co = e α ρ co /(3m 2 p H 2 ) and Ω c = ρ c /(3m 2 p H 2 ), we have 1 = Ω V + Ω co + Ω c . Also from Eq. (5.4) it follows that, if the asymptotic adiabatic regime has been reached, Ω co = γΩ V , and we we obtain Since Ω c ∼ 0.3 today, and γ ≪ 1 in the strong coupling limit we are considering, we must have Ω co ≪ 1 today. The maximum and minimum lengthscales for the instability are still given by Eqs. (5.8) and (5.9), but with ρ replaced by ρ co . Since ρ co is approximately a constant in the strong coupling limit, these lengthscales are also constants. If the parameters of the model are chosen so that Ω c ∼ 1 today, then The evolution equations for the fractional density perturbations δ j = δρ j /ρ j in the adiabatic limit on subhorizon scales are given bÿ where the effective Newton's constants are given by Eq. (C1) with α b set to zero. Writing this out explicitly we obtain The condition for the instability to operate is that the timescale associated with the second term on the right hand side of Eq. (5.21) be short compared with H −1 , or Now the effective mass for this model is given by m 2 eff = β 2 C 2 m −2 p (V + e α ρ co ) = 3β 2 C 2 m −2 p H 2 (Ω V + Ω co ). Substituting this into Eq. (5.22) and using Eq. (5.17) gives the criterion k/(aH) ≫ 1. Therefore the instability should operate whenever modes are inside the horizon and in the range of scales (5.18).
These expectations are confirmed by numerical integrations. In figure 1 we present a numerical analysis of such a two component model. We consider an exponential potential with λ = 2 and strong coupling with C = −20, and typical cosmological parameters are assumed, H 0 = 70 km s −1 Mpc −1 , Ω b = 0.05, Ω c = 0.2, Ω co = 0.05, and Ω V = 0.7. We fix initial conditions of φ/m p = 10 −3 andφ = 0 at a = 10 −10 , although the dynamical attractor renders the final evolution largely insensitive to these choices, we have checked that φ/m p (a = 10 −10 ) = 10 −30 − 1 give the same evolution. We neglect the effect of the coupled CDM component peculiar velocity in the initial conditions (it is many orders of magnitude smaller than the density perturbation) and assume that the coupled and uncoupled CDM components have the same initial fractional density perturbations δ c = δ co , fixed by the usual adiabatic initial conditions. As shown in figure 1, the background evolution is entirely consistent with a ΛCDM like scenario, with w eff = −0.69 today and an asymptotic equation of state at late times, given by (5.11), w eff = −1/(1+γ) = −0.89. The large coupling drives the evolution to an adiabatic regime at late times, with an adiabatic sound speed c 2 a → −1/(1 + γ) as in (5.6). This drives a rapid growth in over-densities once in the adiabatic regime so that although consistent with structure observations at early times, they are inconsistent once the accelerative regime has begun.
In summary, these models provide a class of theories for which the background cosmology is compatible with observations, but which are ruled out by the adiabatic instability of the perturbations.

C. Chameleon models
Next we study the so-called chameleon models [19,20] defined by the potential where M is a mass scale and n > 0 and λ are dimensionless constants, together with the coupling function (5.2). In these models it has been previously shown that the adiabatic regime is achieved in static solutions describing macroscopic bodies like the Earth, and that cosmological solutions in the adiabatic regime provide good models of dark energy [37,40,41]. We now study under what conditions these models are unstable. The effective potential is, from Eq. (2.11), 24) and solving for the local minimum of this potential yields the relation between φ and ρ in the adiabatic regime: x n+1 e x = ρ crit ρ , (5.25) Here x = −βCφ c (ρ)/m p is dimensionless and the critical density is As before the existence of a local minimum in the effective potential requires C to be negative. We shall restrict attention to the regime ρ ≫ ρ crit (5.27) since for models of dark energy ρ crit will be of order the present day cosmological density. In this regime the solution to Eq. (5.25) is approximately (5.28) The corresponding effective mass parameter is and since C is negative, we see that this model is always unstable in the adiabatic regime. Inserting the effective mass (5.29) and the sound speed squared (5.31) into Eq. (4.12) we obtain the range of spatial scales L over which the instability operates for a given density ρ: where We see that the range of unstable lengthscales is nonempty only if which as before is equivalent to the strong coupling condition (4.15). Note that the first of the two inequalities in Eq. (5.32) is equivalent to the "thin-shell condition" of Ref. [20] when the background value of the scalar field can be neglected. As for the exponential models of Sec. V A, the effective equation of state w eff is close to −1 in the adiabatic regime for large coupling, assuming ρ ≫ ρ crit . Therefore these models do not give an acceptable background cosmology for the matter dominated era in their adiabatic regime. However, one can construct two component models analogous to those in Sec. V B using the chameleon potential. Those models give an acceptable background cosmology, but are then ruled out by the adiabatic instability.

D. Mass Varying Neutrino ("MaVaN") models
The impact of adiabatic instabilities has been discussed extensively in the context of MaVaN models, in which the light mass of the neutrino and the recent accelerative era are twinned together through a scalar field coupling [22, 23,24]. The adiabatic instability was shown to be a concern in these models with the implication of forming compact localized regions of neutrinos after undergoing dramatic adiabatic collapse [15].
The action for the MaVaN models is of the form (2.2) with p = 0. Therefore the analysis in the earlier sections of this paper, which assumed p = 2, does not directly apply to these models. However it is straightforward to generalize our analysis to cover this case. Consistent with the general discussion in sections III-IV, the instability is present in these theories if the coupling is strong. Here we discuss examples of models which evade or are subject to the instability on cosmological scales.
Recently, Ref. [17] discussed a MaVaN scenario with a logarithmic potential and scalar dependent mass, where ξ is a constant and m ν0 and φ * are the current values of the neutrino mass and the scalar field respectively. It was found in [17] that this model can exhibit an instability in growth for an otherwise cosmologically viable background solution. We find, however, that this model can also allow stable solutions for identical fractional densities today as those studied in [17], for a wide range of parameter values. If the coupling, m ′ ν /m ν , is not large compared to m −1 p at late times, the evolution never enters the adiabatic regime on cosmic scales. This translates in this model to φ * not being significantly less than m p .
The MaVaN model evolves according to a coupled Klein Gordon equation, ) . (5.37) Assuming no chemical potential, and g spin states per neutrino species with momentum p and mass m ,the neutrino density and pressure are given by with q ≡ ap/k B T 0 ν andm ν ≡ m ν c 2 /k B T 0 ν . In the relativistic regime, withm ≪ 1, the potential is negligible and the driving term, on the right hand side of (5.37), can be calculated by doing a Taylor expansion to first order inm ν , so that where ρ 0 ≡ 7π 2 g(k B T 0 ν ) 4 /240 would be the relativistic neutrino energy density per neutrino species today with temperature T 0 ν . Putting (5.43) into (5.37) and neglecting the potential we find a power law attractor φ ∝ τ x with x = 0.5. The normalization of φ is wholly specified in the attractor by (5.37). Writing φ = φ i (τ /τ i ) 0.5 , and a ∝ τ p , with p = 2/(1 + 3w eff ), we find attractor while the neutrino is relativistic allowing late time evolution to be independent of initial conditions. The vacuum expectation value of the scalar field, if the field becomes adiabatic, (dashed lines) is also shown. For φ * 10 −2 mp the scalar field does not enter an adiabatic era on cosmological scales before now, and growth of perturbations remains wellbehaved. For smaller φ * , for example φ * ∼ 10 −3 mp shown, the evolution is adiabatic at late times, similar to that discussed in [17]. When the neutrino is non-relativistic, if we again neglect the potential, The Klein-Gordon equation has a solution φ ∝ τ x with x = (2 − p)/3, tending towards a cessation of growth in φ in the matter dominated era. If the neutrinos are non-relativistic (but, of course, as-suming that they decoupled when they were relativistic) At late times, the complete effective potential is relevant, for which there exists a minimum at positive φ given by We modified CAMB [39] to investigate the evolution numerically. In figure 2 we show the scalar field and neutrino mass evolution for ξ = 10 20 m −1 p , φ * = 1.8m p , and m ν0 = 0.312eV (giving Ω ν = 0.02), Ω b = 0.05, Ω c = 0.23, H 0 = 70kms −1 Mpc −1 . In the bottom left hand figure, the numerical evolution of scalar field is shown for φ * = 1.8m p along with two smaller values φ * ∼ 10 −3 m p and 0.3m p for which the minimum of the effective potential is steeper. The background evolution obeys the φ ∼ τ 0.5 attractor and normalization in (5.44), rendering it largely independent of the initial conditions. When the neutrino becomes non-relativistic the evolution slows and for small values of φ * starts to track the VEV, as discussed in [17]. However for larger values of φ * , the VEV is not reached until later times, and in addition, the effective potential minimum is shallow enough that the scalar field does not get fixed at the VEV automatically, enabling well-behaved growth. The bottom right hand figure shows the resultant matter power spectrum, which is very similar to a fiducial ΛCDM model with the same Hubble factor, baryon and neutrino density today when normalized to the same amplitude at large scales.

VI. CONCLUSIONS
If dark energy and dark matter are to fit into a coherent fundamental physics framework then there are likely couplings between them. Such couplings may have far reaching macroscopic implications on scales ranging from the solar system up to cosmological horizon. As we have discussed, these effects may lead to observationally distinctive characteristics, which may allow us to tease out the nature of the dark sector. However they may also give rise to catastrophic instabilities with which we may constrain the class of physically viable dark energy models.
In this paper we have considered such theories in which there exists a nontrivial coupling between the dark matter sector and the sector responsible for the acceleration of the universe.
We have comprehensively analyzed an instabilitycharacterized by a negative sound speed squared of an effective coupled dark matter/dark energy fluid -that exists whenever such theories enter an adiabatic regime in which the scalar field faithfully tracks the minimum of the effective potential, and the coupling strength is strong compared to gravitational strength. The adiabatic regime occurs when the relaxation time scale associated with the scalar field is much shorter than the Hubble time. We have demonstrated how this instability can be viewed from the alternative perspectives of the kinetic theory of dark matter and as a Jeans instability associated with modified Newton's constants.
We have established the conditions under which the adiabatic instability occurs, finding a condition on the coupling, |α ′ |m p ≫ 1, and have identified the time and length scales over which the instability is active in a given setting governed by the matter density. These length scales can differ greatly dependent on whether one is considering galactic or cosmic densities.
Our work builds on previous analyses of MaVaN [15,17,22,23,24], chameleon [19,20], and general coupled scenarios [16]. Our numerical analyses of coupled CDM and MaVaN models reinforce our analytic findings. We show that the regime of adiabatic behavior agrees with the predictions based on coupling strength and length scale conditions mentioned above. In the appropriate limits, our results reduce to the previous findings in a number of cases, while in several other cases we have provided corrected results. In particular, we have shown that stable MaVaN models exist that evade the instability, if m ′ ν /m ν is not large compared to m p .
In the body of this paper we have neglected the effect of baryons and the other matter species, assuming that their density is much smaller than that of the dark matter. A natural question is does the instability persist in regimes where the density of visible matter density is comparable to or larger than that of the dark matter -one might expect the instability to be killed by the pressure of the visible matter. In general the answer depends on the coupling function α b in the action (2.1), which governs the coupling of the scalar field φ to normal matter. The instability persists in two specific cases, (i) α b = 0 where normal matter is uncoupled from φ, and (ii) α b = α c , the case usually considered in which dark matter and normal matter couple the same way to the scalar field, respecting the equivalence principle.
To derive this result, we generalize the derivation given in the body of the paper from a single, pressureless fluid to a set of N non-interacting fluids, each of which can have a pressure. We start from the action (2.1), and we define the Jordan frame metric that the jth sector couples minimally to byḡ j ab = e 2αj(φ) g ab . (A1) The stress energy tensorT ab j of the jth sector is defined by and we assume it has the form of a perfect fluid Hereū a j is the 4-velocity which is normalized according tō g j abū a jū b j = −1. This perfect fluid assumption requires that the matter in the jth sector be barotropic, ie that its pressure be determined uniquely by its density, which in many regimes is a good approximation.
The equation of motion of the jth fluid is where∇ j a is the derivative operator determined by the metricḡ j ab . The components of this equation perpendicular to and parallel to the 4-velocity arē whereā a j =ū b j∇ j bū a j is the Jordan frame 4-acceleration, and∇ j a (ρ j +p j )ū a j =ū a j∇ j apj .
Next, we conformally transform these fluid equations using Eq. (A1) to write them in terms of the Einstein frame metric g ab and the Einstein-frame normalized 4-velocities u a j = e αjū a j which satisfy g ab u a j u b j = −1. This gives where a a j = u b j ∇ b u a j and Next, we consider the non-relativistic limit of the fluid equations (A7) and (A8) together with the adiabatic limit of the scalar field equation (2.4). We also neglect selfgravity in the Einstein frame, taking g ab ≈ η ab . This approximation should be valid on sufficiently small spatial scales. We get where α j = α j (φ) and φ = φ(ρ k ) is given by the equation Next, we switch to using the rescaled density variables We linearize the resulting equations about the back- In other words we assume that the dark matter and visible matter densities are constant and that the fluids are not in relative motion. We define α j0 = α j (φ 0 ) and α ′ j0 = α ′ j (φ 0 ), and look for a solution From Eq. (A11) we obtain that δφ = j χ j δρ j , where Next we assume that all the variables are proportional to exp [ik · x − iωt], and that the velocities are proportional to k. This assumption will yield one particular set of modes of oscillation of the coupled fluids, but these modes will contain the instability if it is present. The resulting eigenvalue equation for ω 2 is where the matrix Γ jk , which plays the role of the effective squared sound speed, is is the physical (Jordan frame) squared sound speed of the jth fluid. From the eigenvalue equation (A16), it follows that the system will be stable if and only if all of the eigenvalues of the matrix Γ jk are real and positive. For the case of a single fluid withc 2 j = 0, this criterion reduces to the criterion (4.3) derived in the body of the paper. In general, there is a competition between the positive squared sound speeds of the fluids in the first term in Eq. (A17), and the (possibly) negative squared sound speeds coming from the interaction with the scalar field given in the second term.
We now specialize to two fluids, a dark matter fluid with zero sound speed, and a fluid describing the visible matter. Taking j = 1 = c for the dark matter (CDM) and j = 2 = b for the visible matter (baryons), and defining gives for the matrix (A20) Now if the system is stable, then both of the eigenvalues λ 1 and λ 2 of the matrix must be real and nonnegative, and hence both the determinant λ 1 λ 2 and the trace λ 1 + λ 2 of the matrix must be nonnegative. Conversely, if either the determinant or the trace of the matrix is negative, then the system is unstable. The determinant is where and where we have used Eqs. (A15) and (A19). In Eqs.
(A21) and (A22), φ can be taken to be the function of ρ c and ρ b given by [cf. Eq. (A8)] Equations (A21) -(A23) allow us to determine which values of ρ c and ρ b satisfy the sufficient condition det Γ < 0 for instability, given the functions V (φ), α c (φ) and α b (φ). Note that α c (φ) was denoted by α(φ) in the body of the paper. We now consider two special cases. If α ′ b = 0, so that the visible matter is not coupled to the scalar field, then eigenvalues of the matrix Γ are justc 2 b and ν c χ c . The second of these eigenvalues coincides with the effective squared sound speed computed above in Eq. (4.2). Thus we recover the results of the body of the paper.
The second special case is when α b = α c . In this case the instability criterion det Γ < 0 again reduces to the criterion (4.2) computed earlier for a single fluid, with the modification that the effective sound speed squared is now a function of the total density ρ = ρ c +ρ b rather than just of ρ c . If this total density is in the unstable regime, then the instability will persist despite the pressure of the visible matter.
Finally, one can also consider the effect of the scalar field on just the baryons, neglecting the dark matter. The sound speed of the baryons gets an additional term, but it is a small correction unless the coupling is large, α ′ b m p ≫ 1, and observational tests of general relativity in the Solar System require α ′ b m p ≪ 1.

APPENDIX B: KINETIC THEORY TREATMENT OF INSTABILITY
In this appendix we describe dark matter using the collisionless Boltzmann equation, and specialize to the nonrelativistic regime and to lengthscales small enough that self-gravity can be neglected. We show, first, that the instability is generic, occurring for any velocity distribution function with sufficiently small velocity dispersion, and, second, that for Maxwellian distributions the fluid is stabilized by free streaming once the velocity dispersion becomes larger than a critical value.
Starting from the general action (2.1), we specialize to the non-relativistic limit neglecting self-gravity (g ab ≈ η ab ), which will be valid on sufficiently small spatial scales. The Jordan-frame metric is then We assume that dark matter is composed of noninteracting, non-relativistic particles of mass µ. We denote by f (t, x, v) the one-particle distribution function, normalized so that the number of particles in the volume element d 3 x and in the velocity region d 3 v is Now the physical (Jordan-frame) volume element is e 3α d 3 x, so the Jordan frame mass densityρ [cf. Eq. (2.3) above] is given byρ = µe −3α d 3 vf. From Eq. (2.7) the rescaled density variable ρ is then given by (B1) With this notation the collisionless Boltzmann equation is [cf. Eq. (A10) above withp j = 0] In the adiabatic limit, the variable α in this equation is given by α = α(φ), where φ is given in terms of ρ by Eq. (3.4), and ρ is given in turn in terms of f by Eq. (B1). We now linearize the Boltzmann equation about a homogeneous background solution, taking f = f 0 (v) + δf (t, x, v). and taking δf ∝ exp[−iωt + ik · x]. Using the identity dα/dρ = c 2 s /ρ from Eq. (4.2), and integrating over v, gives the dispersion relation Here ρ 0 is the background density. Without loss of generality we take k to be in the z-direction, assuming f 0 is isotropic, and we define This distribution function is normalized so that dvF = 1, and the dispersion relation can now be written as Given the distribution function F (v) we can solve this equation to obtain ω = ω(k), which is complex in general. The integrand in Eq. (B5) contains a singularity at v = ω/k. The derivation given here is incomplete since it does not provide a specification for how to deal with the singularity. However, just as in plasma physics [42], it is possible to give an alternative derivation based on Laplace transforms. That alternative derivation yields the following specification: the integral over v must be taken over the Landau contour, which runs along the real axis if Im(ω) > 0, but dips below the real axis to encircle the pole at v = ω/k if Im(ω) ≤ 0.
We now specialize to the adiabatic regime where c 2 s < 0, and we write c 2 s = −β 2 . We write the complex frequency in terms of its real and imaginary parts, ω = ω r + iω i , and we look for a solution ω = ω(k) of the dispersion relation (B5) with ω i > 0, corresponding to an unstable mode. For such a solution the Landau contour is along the real axis which simplifies the analysis.
If the velocity dispersion of the distribution F (v) is small, we can expand the denominator of the integrand in Eq. (B5) as a power series in vk/ω. Integrating by parts and solving the resulting equation for ω 2 gives where σ 2 = dvv 2 F (v). This equation has a solution with positive imaginary part, ω i > 0, consistent with the assumption used in its derivation. Therefore the instability is generic, present for any velocity distribution, as long as σ ≪ β. Equation (B6) also suggests that the instability will be removed when σ gets large, as argued in Ref. [15], since then the positive second term in the square brackets will overcome the negative first term. However, Eq. (B6) is only valid in the regime σ ≪ β, and so to investigate stabilization we must use an alternative method of computation.
Returning to the general dispersion relation (B5), we obtain from its real and imaginary parts the equations and We now specialize to the Maxwellian velocity distribution For this case Eq. (B8) has a unique solution for ω r , namely ω r = 0. This can be seen by splitting the integral into v > 0 and v < 0 contributions and substituting v → −v in the v < 0 term. This yields and the result follows since the integrand is everywhere positive. Equation (B7) now simplifies to where we have definedv = v/σ and κ = ω i /(kσ). Equation (B11) determines κ and thence ω i as a function of σ/β. When κ is large, expanding the factor in round brackets in the integrand as a power series inv/κ gives the result (B6) above. As κ decreases, σ/β increases, until as κ → 0, σ/β approaches the limiting value σ/β = 1. It can be seen that there are no solutions to Eq. (B11) with σ/β > 1, since the right hand side is a monotonic function of κ 2 . Therefore, whenever 15 there are no unstable modes with ω i > 0. In other words, the instability has been removed by the damping process associated with the finite dispersion (Landau damping, also called free streaming). The perturbation to the distribution function is proportional to the integrand in Eq. (B11). Therefore the nature of the damping is that interaction with the growing mode moves some particles from velocities v ∼ β to larger velocities, removing energy from the mode.

APPENDIX C: EFFECTIVE NEWTON'S CONSTANT
In this appendix we derive the formula (C1) for the effective Newton's constant for the theory (2.1).
Since Newton's constant G is dimensionful, we need to define the system of units used in measurements of G, as our final result will depend on the system of units chosen. Here we are not concerned with changes of units in the usual sense where the ratios between the old and new standards of mass, length and time are constants, independent of space and time. Rather, we are concerned with changes in the operational procedures for how units are defined, for which the ratios between the old and new standards can vary with space and time, as discussed by Dicke [43]. Changes of units of this type encompass 15 Note that the critical value of velocity dispersion σ = β is a factor of √ 3 larger than indicated by the approximate formula (B6) which was used in Ref. [15]. conformal transformations of the metric, but can be more general. The necessity of specifying a system of units before discussing the space-time variation of dimensionful constants of nature is discussed in detail in Duff [44].
We choose to use systems of units which are determined completely by physics in the visible sector (which we label by j = b for baryonic), and are not determined by gravitational physics or by physics of the dark matter sector. For example, this is true of SI units. Alternatively one could take one's standards of mass, length and time to be the mass, size, and light travel time across a Hydrogen atom. Any system of units of this type will give rise to the formula (C1) below for G.
We note that alternative choices are in principle possible. For example, one could have a small fiducial black hole, and one could transport it adiabatically around the Universe to act as a standard. Then one could define the standard of mass to be the mass of the black hole, while defining the standards of length and time in terms of Hydrogen atoms as above. In this system of units G would vary in space and time, but in a manner different to that described by the formula (C1) below, and Planck's constant would also vary in space and time. Another possibility would be to define the units of mass in terms of the mass of the dark matter particle; this would also yield a different formula for G whenever the dark matter coupling function α c (φ) differs from that of the baryons α b (φ).
We next discuss the quantities on which the formula for G depends. First, for the theory (2.1) with several different matter sectors, Newton's constant becomes a matrix whose i,j element governs the strength of the gravitational interaction between sector i and sector j. Second, G will depend on the scalar field φ, and through this dependence become a function of space and time. Here φ should be thought of as a background value of the scalar field; gravitational interactions that act as perturbations to this background, for which the perturbation to the scalar field can be treated linearly, have a strength described by the formula (C1). Finally, since the scalar contribution to the gravitational force will in general have a Yukawa profile, G will also depend on a spatial wave vector k, which is measured in the units discussed above.
Our formula for Newton's constant is Here on the right hand side, G is a constant, related to m 2 p by G = 1/(8πm 2 p ), and the effective mass m eff (φ, ρ) is defined by where the effective potential is given by Eq. (2.11). In the adiabatic regime this effective mass reduces to the mass (3.5). Several different aspects of this formula have appeared before in the literature. The overall prefactor of e 2α b is well known, and the formula for the case of one matter sector has been previously derived in the context of cosmological perturbation theory [45]. The 1 in the square brackets describes the exchange of a graviton, and the second term in the square brackets describes the exchange of a scalar quantum, which couples to particles in the sector i with an amplitude proportional to α ′ i (φ). The scalar coupling term vanishes in the long wavelength limit k → 0 if the scalar field has a finite effective mass m eff .
We now turn to the derivation of the formula (C1). For the derivation it will be convenient to use a different system of units, which are defined as follows. The standards of mass, length and time are defined to be the mass, size and light travel time across a small fiducial black hole, which is transported adiabatically around the Universe to act as a reference. Equivalent units 16 can be defined by demanding that the speed of light c and Planck's constant be unity, and that lengths and times are measured using the Einstein frame metric g ab . At the end of the derivation we will transform back to the visible-sectorbased units.
Starting with the action (2.1), we perform the following steps. We specialize the matter action in the jth sector to be that of a collection of point particles with masses m (C3) We specialize to the Newtonian limit, so that the Einstein-frame metric g ab is of the form −(1 + 2Φ)dt 2 + (1 − 2Φ)dx 2 . We write the scalar field as φ + δφ, and expand to quadratic order in δφ and Φ. We assume that there is a background mass density ρ, so that the effective mass of the scalar field is given by Eq. (C2). Finally, we compute from the action the energy E for a static configuration, up to an additive constant. This yields where