# Dilute and dense axion stars

###### Abstract

Axion stars are hypothetical objects formed of axions, obtained as localized and coherently oscillating solutions to their classical equation of motion. Depending on the value of the field amplitude at the core , the equilibrium of the system arises from the balance of the kinetic pressure and either self-gravity or axion self-interactions. Starting from a general relativistic framework, we obtain the set of equations describing the configuration of the axion star, which we solve as a function of . For small , we reproduce results previously obtained in the literature, and we provide arguments for the stability of such configurations in terms of first principles. We compare qualitative analytical results with a numerical calculation. For large amplitudes , the axion field probes the full non-harmonic QCD chiral potential and the axion star enters the dense branch. Our numerical solutions show that in this latter regime the axions are relativistic, and that one should not use a single frequency approximation, as previously applied in the literature. We employ a multi-harmonic expansion to solve the relativistic equation for the axion field in the star, and demonstrate that higher modes cannot be neglected in the dense regime. We interpret the solutions in the dense regime as pseudo-breathers, and show that the life-time of such configurations is much smaller than any cosmological time scale.

^{†}

^{†}preprint: MCTP-17-20A

^{†}

^{†}preprint: MIT-CTP/4949

^{†}

^{†}preprint: NORDITA-2017-112

## I Introduction

The QCD axion Weinberg (1978); Wilczek (1978); Shifman et al. (1980); Kim (1979); Zhitnitsky (1980); Dine et al. (1981); Preskill et al. (1983); Abbott and Sikivie (1983); Dine and Fischler (1983) arising within the Peccei-Quinn solution of the strong CP-problem Peccei and Quinn (1977a, b) is one of the best motivated dark matter candidates. Other bosonic dark matter candidates include axion-like particles Arias et al. (2012) emerging in many extensions of the Standard Model, especially in string theory compactifications Svrcek and Witten (2006); Arvanitaki et al. (2010); Ringwald (2014); Halverson et al. (2017).

If bosons comprise the dark matter of our Universe, they could form dense (with respect to the average dark matter density) clumps called boson stars. Kaup (1968); Ruffini and Bonazzola (1969), or axion stars in the specific case of axion dark matter. (Here “star” is used to denote an object sustained by hydrostatic equilibrium, whether or not it emits light.) Such objects have been long studied Kaup (1968); Ruffini and Bonazzola (1969); Das (1963); Feinblum and McKinley (1968); Teixeira et al. (1975); Colpi et al. (1986); Seidel and Suen (1991); Tkachev (1991); Kolb and Tkachev (1993, 1994), and recently there has been revived interest Braaten et al. (2016a); Eby et al. (2016a, b); Levkov et al. (2017); Helfer et al. (2016); Braaten et al. (2017, 2016b); Bai et al. (2016); Eby et al. (2017); Desjacques et al. (2017).

In this article, we study the stability of axion stars as a function of the amplitude of the axion field at the core of the star . Our results apply to the full range of axion masses for which QCD axions can comprise all of the dark matter. We identify three distinct branches of axion stars, distinguished by the field amplitude at the core, which in turn determines the density of the star. We should keep in mind, that the axion is a periodic field with amplitude effectively restricted to the domain .

For small field values with the axion mass, the axion field only probes the harmonic part of the potential, and it can be treated as a free field. In this regime, self-gravity is balanced by the kinetic pressure arising from the uncertainty principle. We call this the dilute axion star branch. We reproduce the previous findings in the literature for the mass-radius relationship, , where and are the radius and mass of the star, respectively. In this regime, the configuration is stable against perturbation: For a given mass , stars are pulled back to the equilibrium radius if they expand because then the (attractive) self-gravity is stronger than the (repulsive) kinetic pressure; conversely, if they are perturbed to smaller radii, they expand because kinetic pressure becomes stronger than self-gravity.

For configurations with , self-interactions cannot be neglected anymore, although the amplitude is still comparatively small, . For QCD axions, the lowest order self-interaction is an attractive quartic term. For amplitudes , the attractive quartic self-interaction is stronger than gravity, which is negligible in this regime. In this critical branch, we find solutions when the quartic self-interaction balance the kinetic pressure with mass-radius relation . Note, that this relation implies that axion stars become lighter with growing density, such that they always have masses in this branch. However, the solutions are unstable against perturbations: for a given mass , stars expand when perturbed to radii larger than the equilibrium value since the quartic self-interactions are weaker than the repulsive pressure. Eventually the configuration relaxes to the dilute, interaction-free regime described in the previous paragraph. Conversely, if configurations are perturbed to radii smaller than the equilibrium value, the quartic interaction is too strong to be balanced by the pressure and the star collapses to even higher densities.

It has recently been pointed out, that new stable configurations, called dense axion stars, are obtained when the amplitude of the axion field in the core reaches Braaten et al. (2016a). For such amplitudes, the axion field scans the full non-perturbative axion potential, and self-interactions must be taken into account to all orders. Using the assumption that the axion field in the star is coherently oscillating at a single frequency, as commonly used in the literature, we obtain the mass-radius relation , in agreement with Ref. Braaten et al. (2016a). However, we find that the single-harmonic approximation, which holds in the branches described above, is not accurate for the dense branch. Using a multi-harmonic expansion, we find that higher harmonics are generated with amplitudes comparable to the fundamental mode’s amplitude. Heuristically, the presence of higher harmonics corresponds to the generation of (relativistic) axions by coalescence processes . We find that configurations on the dense branch decay via emission of relativistic axions, with lifetimes of order , which are much shorter than any cosmological timescale.

When , axions stars are short lived solutions of the relativistic equation, elsewhere known as oscillons Dashen et al. (1975); Kudryavtsev (1975); Makhankov (1978); Geicke (1984); Zakharov and Kuznetsov (1986); Gleiser (1994); Copeland et al. (1995); Salmi and Hindmarsh (2012); Mukaida et al. (2017). In the literature, similar objects have also been called pseudo-breathers Hormuzdiar and Hsu (1999), axitons Kolb and Tkachev (1994), or oscillatons when driven by gravity Seidel and Suen (1991); Matos and Ureña-López (2000); Urena-Lopez (2002). Since gravity is negligible in the dense branch, the axion field is described by the Klein-Gordon equation with the QCD chiral potential (the -Gordon equation). There is a large but scattered literature on finding solutions to related equations. For example, in one dimension, assuming a cosine potential leads to the Sine-Gordon equation, which admits localized breather solutions that are not harmonic Ablowitz et al. (1973), i.e. which feature an infinite collection of higher harmonics. In three dimensions oscillons closely resemble the breather solutions of the one dimensional Sine-Gordon equation, but they differ in that they radiate energy and thus decay in a finite lifetime, though slowly relative to the “natural” timescale set by the inverse mass of the particles.

Justifying and expanding upon this concise summary, the remainder of this paper is as follows. In Sec. II we set out the basic equations. In Sec. III we find numerically stable solutions and provide quantitative results for the dilute and the critical axion star branches. In Sec. IV we discuss the dense branch, and analyze the equilibrium and metastability of dense configurations in a relativistic framework. In Sec. V we use the mass-radius diagram to sketch a qualitative storyline for axion stars, and in Sec. VI we summarize and conclude.

## Ii Axion stars

### ii.1 Axion Lagrangian

The axion results from promoting the flavor-neutral CP violating angle of the standard model, , to a dynamical field Wilczek (1978); Weinberg (1978) in the Peccei-Quinn mechanism Peccei and Quinn (1977a, b). The canonical normalization of the dynamical angle requires a new energy scale , the axion decay constant, to define the axion field . In the following we will refer to both and as the axion field. The dynamics of the axion field under the influence of gravity are described by the action

(1) | |||||

where the metric is determined by the Einstein equation for the energy momentum tensor of the axion field . We adopt the axion potential Di Vecchia and Veneziano (1980); Grilli di Cortona et al. (2016),

(2) |

where is the topological susceptibility Grilli di Cortona et al. (2016); Petreczky et al. (2016); Borsanyi et al. (2016) and with the ratio of the up and down quark masses . Note, that the minimum of the potential is at and the maximum at . The axion mass and the quartic coupling constant are defined through

(3) | |||||

(4) |

Assuming spherical symmetry and expanding the metric to linear order about flat space yields the line element

(5) |

where is the gravitational potential, which satisfies the Poisson equation with energy density , and is the differential solid angle. In the following, we rescale time and radius as and , respectively, so that the Lagrangian in Eq. (1) reads

(6) |

where a dot indicates a derivative with respect to the rescaled time, a prime indicates a derivative with respect to the rescaled radius, and . Coupling the Poisson equation with the equation of motion obtained from the Lagrangian density gives

(7) | |||

(8) | |||

(9) |

where with the Planck mass GeV. The energy density is dimensionless, and reaches when and the axion potential saturates. In Eq. (9), we denote the contributions to the energy density from the kinetic, gradient, and potential components separately. Note, that the gradient energy is due to the momentum of the axion arising from the uncertainty principle. So far, the only approximation used is that gravity is weak, .

We anticipate one of the results of this paper, namely that the system can be studied in two different regimes depending on whether the axion field is (the “dilute” and the “critical” axion star regimes) or (the “dense” axion star regime). In the dilute and critical regimes, the axions comprising the star are non-relativistic and the tools described in Sec. II.2 below apply. When , a full relativistic description is needed, as we sketch in Sec. III.3.

### ii.2 Non-relativistic (single harmonic) limit

When the non-relativistic limit applies, the axion mass is the largest energy scale in the problem, so that axion stars oscillate at a frequency very close to the axion mass . Despite non-linear interactions arising from a cosine or a chiral potential precluding axion stars solutions from having one single frequency, for small field configurations , the one-frequency approximation

(10) |

suffices. Here, is the total energy of a constituent axion, in units of the axion mass. We write , where accounts for the contribution from the binding, kinetic and self-interaction energies, while the one accounts for the rest mass energy. In the non-relativistic approximation, we have and .

We further assume that gravity is a weak effect, so that we can drop all terms containing in Eq. 7, except for the term which is of the same order as . We split the potential into a mass term and the self interaction as

(11) |

Inserting the representation in Eq. (10) into Eqs. (7)-(9) and averaging over the period , we obtain

(12) | |||||

(13) | |||||

(14) |

In the last expressions, we have introduced the energy density terms

(15) |

and we have defined the effective self-interaction potential and its first derivative through

(16) | |||||

(17) |

For , Eq. (12) is a Schrödinger equation for the radial eigenfunction with eigen-energy , while the energy density reduces to since the contributions from the gradient term and self-interactions are negligible.

We stress that our procedure, which involves the average over of the equation of motion leads to the same results as what was obtained in Ref. Guth et al. (2015), where the authors neglect the rapidly oscillating terms proportional to powers of . As long as gravity is negligible and the single-harmonic approximation in Eq. (10) holds, Eqs. (12)–(14) are valid even for relativistic axions. We anticipate, that for (most of) the dense branch, gravity is indeed negligible but the single harmonic approximation no longer holds.

### ii.3 Axion potential

We expand the expression in Eq. (2) as

(18) | |||||

(19) | |||||

(20) |

In our numerical calculation, we truncate the sum in Eq. (18) to the first five terms . This attains a precision below 1% with respect to the chiral potential in Eq. (2); this precision is better than the accuracy of the chiral perturbation theory itself. We slightly modify the coefficients so that the truncated potential shows: I) the same minimum , II) the same mass , and III) the same quartic coupling as the full chiral potential in Eq. (2), where the (negative) quantity is related to the axion quartic self-interaction constant as . The numerical values of the corresponding corrected coefficients are given in Table 1 for .

The effective non-relativistic potential in Eq. (16) is

(21) |

where is the Bessel function of the first kind of order zero for the argument . Notice that the cosine potential is recovered in the limit , equivalent to setting , , and all other equal to zero in Eq. (18). The set of Eqs. (12)-(14) has been extensively applied to self-gravitating systems made of bosons. For the case of axions, the free case has been studied in Refs. Seidel and Suen (1991); Widrow and Kaiser (1993); Marsh and Pop (2015) following the seminal work in Refs. Kaup (1968); Ruffini and Bonazzola (1969). The potential expanded to the quartic interactions has been studied in Refs. Chavanis (2011); Chavanis and Delfini (2011); Guth et al. (2015); Cotner (2016). Ref. Braaten et al. (2016a) considers the set of Eqs. (12)-(14) with the cosine potential, using the expression for the energy density (in our notation) , instead of our Eq. (15) obtained from the full energy-momentum tensor. This implicitly neglects contributions from self-interaction and kinetic energy to the energy density, which sources the gravitational potential. As we show below, those contributions to the energy density affect the results for the “dense” branch.

## Iii Numerical results in the single harmonic approximation

### iii.1 Axion star branches

We numerically solve for the radial profile appearing in the set of Eqs (12)-(14), as a function of the frequency . We impose the boundary conditions

(22) |

where is the rescaled energy density at and the core amplitude is the amplitude of the axion field at . We obtain a radial profile via a shooting method, that is by varying the value of the core amplitude until we find a profile that decays as at a sufficiently large . The solution we seek shows no nodes, and corresponds to the lowest energy state for a given value of . See Ref. Davidson and Schwetz (2016) for excited states of an axion star with a quartic potential. We find solutions for all values of within the range (0,1), although the numerics are particularly tricky as we approach . For each value of , we obtain a unique value of the core amplitude and a unique profile. Given the radial profile, we obtain the total mass and the radius of the axion star, the latter defined as the radius containing 90 % of the energy Ruffini and Bonazzola (1969). In Fig. 1, we show the mass-radius relation for three values of GeV.^{1}^{1}1Note, that for GeV some fine-tuning of the misalignment angle is required to avoid overclosure of the Universe Hertzberg et al. (2008); Visinelli and Gondolo (2009). Each point on the line is characterized by a fixed value of and the core amplitude . For increasing value of , we identify three different regimes: the dilute branch (), the unstable critical configurations (), and the dense branch (). For the critical line (red dashed line) and (most of) the dense branch (dashed black line), gravity is negligible. Then we find universal solutions when expressed in terms of the natural units of star mass, , and radius, . However, gravity is relevant in the dilute branch, where solutions depend on the value of through .

### iii.2 Non-relativistic solutions

In this section, we present heuristic arguments explaining the numerical results obtained in the previous Sec. III.1 for the dilute and critical branches where ; see also Chavanis (2011); Chavanis and Delfini (2011) for a similar approach. These branches can be understood in terms of the different contributions to the axion star energy : the gravitational binding energy, the gradient energy, and the (quartic) self-interaction contribution,

(23) |

Here, and are dimensionless parameters which we insert to match the analytical results derived from Eq. (23) with the numerical solution. Estimating the mass of the axion star as

(24) |

we can express the central amplitude as , and the total energy can be rewritten as

(25) | |||||

In the last equality, we have used the scaling property of the Schrödinger-Poisson equation, writing the mass and the radius of the star in terms of dimensionless quantities, and . The natural scale for the mass and the radius of the axion star are then

(26) | |||||

(27) |

where and are respectively the mass and the radius of the Sun. The equilibrium configurations of the axion star can be qualitatively obtained by minimizing the energy density in Eq. (25) with respect to , while fixing the axion star mass or, equivalently, the total number of axions . This gives a quadratic equation whose solutions correspond to the radius of the star for either the dilute branch () or the critical branch (), namely

(28) |

The stability of the solution is determined by the sign of . Solutions in the dilute branch () are stable, while those in the critical branch () are unstable. Matching onto our numerical results from section III.1, we obtain

(29) |

independent of the value of .

The dilute branch of the axion star corresponds to the equilibrium between the gradient energy and gravity. Depending on the value of the decay constant, equilibrium configurations of this type populate the line with negative slope in Fig. 1 with GeV (blue), GeV (green), or GeV (orange), with the mass-radius relation

(30) |

For configurations lying above this equilibrium line, the gravitational pull overcomes gradient pressure, so these configurations contract. On the contrary, configurations lying below the mass-radius line in Eq. (30) are restored to the equilibrium condition by the gradient pressure term. Hence, a restoring force acts to vanish any deviation from the stable equilibrium.

The critical branch, the dashed red line in Fig. 1, corresponds to the balance of the gradient and the quartic self-interaction energy contributions, with mass-radius relation

(31) |

Deviations from this configuration are pushed either further towards the dilute branch or to further contraction and are hence unstable. A solution for the radius of the axion star exists as long as the quantity below the square root in Eq. (28) is positive, that is when the mass of the star is smaller than the critical value

(32) |

which corresponds to the radius and to the core amplitude

(33) | |||||

(34) |

The values of and define the turning point in the top right corner of Fig. 1, corresponding to the transition from the dilute to the critical branch. In the critical branch, a denser solution corresponds to moving along the red dashed line in Fig. 1 towards the bottom left of the figure, with the star contracting and becoming lighter. Since in this branch the core amplitude increases as , non-perturbative dynamics becomes relevant when , or at a typical mass

(35) | |||||

(36) |

These values of and mark the second turning point in the bottom-left region of Fig. 1. For larger values of the core amplitude, the axion field explores the whole chiral potential and a different treatment is needed.

### iii.3 Non-perturbative solution

The axion star solutions found for correspond to a clump of axions whose total mass and radius are larger than the critical values in Eqs. (35) and (36). For such configurations, higher order terms in the attractive self-interacting potential cannot be neglected and a new regime is obtained, often referred to as the “dense” axion star regime in the recent literature Braaten et al. (2016a, 2017). We show the numerical results for the mass-radius relation obtained in the dense branch configuration with the solid black line in Fig. 1. Fitting the curve far from the turning point leads to the relation . This regime corresponds to classically stable configuration with an almost constant density in the inner core. For the mass-radius relation, we have obtained the same power-law exponent (1/3) as in Ref. Braaten et al. (2016a), because such dependence follows from the fact that the solution in the dense branch saturates the QCD potential and leads to a constant density of the star.

However, the structure of our solution differs greatly from what was obtained in Ref. Braaten et al. (2016a). We disagree on their interpretation of the equilibrium of the axion star in the dense branch for three main reasons. I) We have included the self-interactions and the gradient energy terms through Eq. (15). These terms cannot be neglected, as we show in Fig. 2. II) In Ref. Braaten et al. (2016a) the set of equations is solved in the Thomas-Fermi approximation, that is neglecting the Laplacian of appearing on the left-hand side of Eq. (12). III) Most importantly, the single-harmonic approximation in Eq. 10 does not hold in the non-perturbative regime.

In Fig. 2, we show the different contributions to the mass of the axion star, , from the various components in Eq. (9), namely , where , as a function of the core amplitude.

In the () regime shown, the star is in the critical (dense) branch. In the critical branch, the kinetic and potential energies both contribute a factor equal to . This result can be interpreted by the fact that the wave function of the coherent axion field undergoes harmonic oscillations, with the energy density equipartitioned between the kinetic and potential terms. However, as we approach the dense regime, the contribution from the gradient term increases, to the extent that for all three components contribute with a similar magnitude. Thus for dense axion stars the energy density must include all energy contributions. Also, the Thomas-Fermi approximation is not justified since the Laplacian term is crucial for solving Eq. (12) in the whole domain shown in Fig. 2 and 3.

In more detail, the structure of a dense axion star looks as follows. The stellar core is composed of relativistic axions since in that region , although self-interactions are not entirely negligible. As we move out of the core, there is an intermediate region where the self-interactions balance the gradient term. Finally, in the outmost part self-interactions are again negligible .

To further illustrate that the axion field is relativistic in the dense regime, in Fig. 2 we show the axion energy per particle (black solid line), which drops to zero for , due to the fact that self-interactions increase with . Then, the non-relativistic condition , which expresses that the typical momentum of the axion is much smaller than its energy, no longer holds. Fig. 3 also shows this conclusion, since the quantity decreases from being much larger than one to a constant value for which the non-relativistic interpretation does no longer hold. The inequality , or , which holds even in the dense branch, is not sufficient to justify a non-relativistic approach.

In addition, our solution shows that gravity is negligible everywhere inside the star. The gravitational energy density at a distance from the center of the star is , where is the mass enclosed within the radius , so we can write

(37) |

where in the last step we have used the parametrization . Hence, gravity can be neglected for . For , gravity can be safely neglected as long as or , which is the range of parameters considered in this work. However, for dense axion star solutions of larger mass, gravity could eventually become important again for . We do not consider this latter possibility here.

As we have previously discussed, the solutions obtained in the dense branch are not self-consistent because the single frequency approximation in Eq. (10) is not justified on the basis of the findings in Fig. 3. When the amplitude of the axion field becomes , the axion fields probes the full chiral potential and all orders of self-interaction become relevant. Then, higher harmonic modes of the axion field whose frequency is a multiple of the fundamental mode are generated with amplitude comparable to that of the fundamental mode. In the next section we therefore start over from Eq. (7) and perform a multi-harmonic expansion.

## Iv Oscillons

### iv.1 Generalities on the relativistic equation

Based on the findings of the previous Section, axions in the dense regime can be studied using a relativistic approach and ignoring gravity. For simplicity, we derive results for the illustrative case of a cosine potential

(38) |

obtained from the chiral potential Eq. (2) for . In that case, the relativistic equation of motion is the Sine-Gordon equation

(39) |

We wish to identify the oscillon solutions of Eq. (39), namely the solutions that are spatially-localized and time-periodic. Such solutions circumvent Derrick’s theorem Derrick (1964), which states that the scalar field Lagrangian in Eq. (1) expressed in flat space-time does not admit time-independent, finite energy solutions because shrinking a non-zero field configuration effectively reduces the total energy of the system Rosen (1968); Lee (1976); Friedberg et al. (1976a, b); Col (1986). Although the ansatz we used previously, Eq. (10), is not a proper solution for the non-time averaged potential, we expect it to be a reasonable approximation at the transition from the non-relativistic to the relativistic domain when .

There is a long history of searching for oscillons of the Sine-Gordon equation, with the most positive outcome being solutions that last oscillations in units of Bogolyubskiǐ and Makhan’kov (1976); Gleiser (1994); Gleiser and Sornborger (2000); Fodor et al. (2006); Salmi and Hindmarsh (2012). The general consensus is that absolutely stable solutions do not exist, although we know of no definite proof. In any case, is much that we can learn about unstable oscillons from the literature.

For axions in particular, Kolb and Tkachev Kolb and Tkachev (1994) discovered the so called “axitons” when studying the cosmological evolution of the axion field in the dark matter context. They followed the evolution of the Sine-Gordon equation in an expanding Universe in which the axion mass strongly depends on the cosmic time, and identified an instability condition that leads to small clumps of the axion field with large values to disappear in bursts of relativistic axions. This instability, which originates from the attractive quartic self-interaction term, is well known in the condensed matter community and has been recently revisited in Ref. Levkov et al. (2017). In that paper, the authors follow the collapse of a dilute axion star with a mass slightly above the critical value . The axion star solution shows a self-similar collapse that ends when the central amplitude saturates the axion potential. Then, the axion field oscillates for a few times, radiating relativistic axions and relaxing to a small amplitude which is nevertheless larger than the starting value. Such instabilities are triggered for a few times until the central amplitude relaxes to the stability region described above. The simulations include gravity, so that the final state can still be a dilute axion star, but the dynamics of the collapse and the radiation of relativistic axions happens at very small radii where gravity is negligible compared with the self-interactions and gradients.

The simulations in Ref. Levkov et al. (2017) are of considerable phenomenological interest, since in principle the collapse of dilute stars is the most natural mechanism to produce dense axion stars. However, one can address the question of dense axion star stability separately from their possible cosmological origin. For such a task we need other means. A promising approach emerged in Ref. Piette and Zakrzewski (1998), where the authors convert the Sine-Gordon equation into a series of equations with different harmonics.

### iv.2 Beyond the 1st harmonic approximation

A general time-periodic solution can be written in terms of an infinite numerable set of harmonics. Thus we can write our oscillon ansatz as as Alfimov et al. (2000)

(40) |

which, once plugged into the Sine-Gordon Eq. (39), yields a set of coupled equations for the different harmonics,

(41) |

Here, we have introduced the notation

(42) |

The set of Eq. (41) is a generalization of Eq. (10) when higher harmonics other than the fundamental mode are considered; when truncating the sum at we obtain the single harmonic approximation Eq. (12) with .

As an example, we consider the case where we also include the first term beyond the single-harmonic approximation besides the fundamental mode . This gives

(43) | |||||

(44) | |||||

(45) |

where we have approximated the computation of the coefficients and by expanding around , with

(46) |

In fact, the solutions found in Sec. III.3 correspond to the zeroth-order approximation of the full non-linear solution, while solving the set of Eqs. (43)–(44) gives the next-to-leading order contribution.

At , solutions must approach zero, with . In this regime, the higher harmonic must satisfy

(47) |

which, in the range , is an oscillatory solution in space with wavelength . Thus, if , the axion star configuration radiates energy away through the third harmonic, with a contribution that increases with the total energy of the axion . Equation (47) can in principle be used to compute the lifetime of the axion star in the dilute branch, where the third harmonics is a very small perturbation of the exact solution. For values of , there is no radiation solution at and higher harmonics have to be taken into account through Eq. (41). In Fig. 4, we solve the set in Eqs. (43)–(44) for a given frequency , with , using a shooting method to obtain the initial conditions for and at that satisfy . We stress that the amplitudes of the 1st and 3rd harmonic are of the same order of magnitude everywhere in the star, demonstrating that the single harmonic approximation sufficient for the case of dilute axions stars does not suffice for the description of the dense regime.

## V Discussion

We have shown that when , axions are relativistic and axion stars enter the dense branch regime where the configuration behaves as a meta-stable oscillon of the -Gordon equation, with a characteristic lifetime. For a free field, bosons stream away from the oscillon core of a star of radius , with a lifetime Gleiser and Sicilia (2009) and with a radiation spectrum peaking at with width . Including a more realistic non-linear self-interaction potential modifies the spectrum by lowering the peak frequency at a lower value , with a new width . Following Ref. Gleiser and Sicilia (2009), an oscillon forms if the two spectra do not significantly overlap, that is when

(48) |

The computation of the oscillon lifetime for a quartic self-interaction has been addressed in Refs. Gleiser and Sicilia (2008, 2009), where the relatively long lifetime (on the scale of the intrinsic timescale ) of oscillons is explained by the relatively small overlap between the oscillation frequencies. Following this method, we estimate of the lifetime of an oscillon for a cosine potential as

(49) |

where we have used the parameters , , and , following Refs. Gleiser and Sicilia (2008, 2009) with the axion Lagrangian in Eq. (1) and a Gaussian ansatz for the radial wave function. In short, the energy of an oscillon is described by its radius and amplitude, and damped oscillations in the oscillon develop along the line of constant minimum energy Gleiser (1994); Copeland et al. (1995). We performed an independent check of these results by using the solutions of the time-independent Eq. (12) in the dense branch as initial conditions which we time-evolve with the Sine-Gordon Eq. (39), as prescribed in Ref. Piette and Zakrzewski (1998). Although this initial wave function is not a proper solution to the Sine-Gordon equation, our numerical solutions yield breather solutions. For a period as considered above, we find that the solution decays after . This result is of the same order of magnitude as what we obtained using Eq. (49)

(50) |

The fact that pseudo-breathers exist has been shown in Ref. Hormuzdiar and Hsu (1999), where the existence of a finite life-time solution to the Sine-Gordon equation has been related to the singular behavior of the solution at zero, when an oscillating function has been imposed as the boundary condition of the solution at infinity. Pseudo-breathers are ultimately decaying states, as discussed in length in Ref. Kolb and Tkachev (1994), where it is found numerically that such solutions are unstable and fragment into smaller clumps. The dynamics and the initial conditions considered in Ref. Kolb and Tkachev (1994) are however different from ours, since the authors consider a cosmological evolution of the axion field with “white noise” initial conditions, and included the Hubble rate in the equation of motion.

Recently, Helfer et al. Helfer et al. (2016) have studied the stability of axion stars including gravity and non-linear effects, finding that stable dense profiles may be possible when , the exact value depending on the axion star mass. In any case, the energy scale involved is well above the scales we consider here. For values of below this critical value, the axion star either collapses to a black hole or dissolves by the emission of relativistic particles, consistently with the puffing out obtained in Ref. Levkov et al. (2017) and in this work.

## Vi Conclusions

In this paper, we have discussed the properties of axion stars for all allowed values of the core amplitude of the axion field . In particular, we have discussed how classically stable solutions can arise from the interplay between self-gravity, axion self-interactions, the pressure due to the Heisenberg uncertainty principle, and the kinetic energy. Using assumptions commonly made in the literature, we have obtained a set of equations describing coherent axion field oscillations inside the axion star in the single-harmonic approximation. For small core amplitudes , we confirmed known results for axion stars in the dilute and critical branches, and provided a heuristic interpretation of those results from first principles.

For , the “dense” regime, we recover similar results to those in Ref. Braaten et al. (2016a) when using the single-harmonic approximation, in particular, the mass radius relation . However, we argue that the single-harmonic approximation does not hold for the dense regime and thus a different approach is needed, taking into account higher harmonics. In the end, we arrive a very different physical interpretation of the dense regime. We find gravity to be negligible for . Dense axion stars should be solutions to the Sine-Gordon (or -Gordon) equation describing the axion field inside the star. We computed the lifetime of dense configurations using both the semi-analytical procedure described in Gleiser (1994); Copeland et al. (1995); Gleiser and Sicilia (2008, 2009) and by using our single-harmonic solutions as initial conditions, which we time-evolved numerically using the Sine-Gordon equation as prescribed in Piette and Zakrzewski (1998); Alfimov et al. (2000). Both methods yield comparable lifetimes of order , much shorter than any cosmological time scale.

We conclude that if dense axion stars can be formed, they would immediately (on cosmological scales) radiate relativistic axions and decay. Since axion stars in the critical branch are unstable against perturbations and either expand to stable dilute configurations or contract to the dense branch and subsequently decay, stable axion stars with mass appear implausible.

## Additional Note

###### Acknowledgements.

We would like to thank Eric Braaten, Katherine Clough, Malcolm Fairbairn, Oleg Gnedin, Thomas Helfer, Z̆elimir Marojević, David J. E. Marsh, and Scott Tremaine for the useful discussions and comments that led to the present work.SB and LV would like to thank the University of Michigan and the Massachusetts Institute of Technology, where part of this work was conducted, for hospitality.

SB, KF, and LV acknowledge support by the Vetenskapsrådet (Swedish Research Council) through contract No. 638-2013-8993 and the Oskar Klein Centre for Cosmoparticle Physics. KF acknowledges support from DoE grant DE-SC007859 and the MCTP at the University of Michigan. JR is supported by the Ramon y Cajal Fellowship 2012-10597, the grant FPA2015-65745-P (MINECO/FEDER), the EU through the ITN “Elusives” H2020-MSCA-ITN-2015/674896 and the Deutsche Forschungsgemeinschaft under grant SFB-1258 as a Mercator Fellow. FW’s work is supported by the U.S. Department of Energy under grant DE-SC0012567, the European Research Council under grant 742104, and the Vetenskapsrådet (Swedish Research Council) under Contract No. 335-2014-7424.

## References

- Weinberg (1978) S. Weinberg, Phys. Rev. Lett. 40, 223 (1978).
- Wilczek (1978) F. Wilczek, Phys. Rev. Lett. 40, 279 (1978).
- Shifman et al. (1980) M. Shifman, A. Vainshtein, and V. Zakharov, Nuclear Physics B 166, 493 (1980).
- Kim (1979) J. E. Kim, Phys. Rev. Lett. 43, 103 (1979).
- Zhitnitsky (1980) A. R. Zhitnitsky, Sov. J. Nucl. Phys. 31, 260 (1980).
- Dine et al. (1981) M. Dine, W. Fischler, and M. Srednicki, Phys. Lett. B104, 199 (1981).
- Preskill et al. (1983) J. Preskill, M. B. Wise, and F. Wilczek, Phys. Lett. B120, 127 (1983).
- Abbott and Sikivie (1983) L. F. Abbott and P. Sikivie, Phys. Lett. B120, 133 (1983).
- Dine and Fischler (1983) M. Dine and W. Fischler, Phys. Lett. B120, 137 (1983).
- Peccei and Quinn (1977a) R. D. Peccei and H. R. Quinn, Phys. Rev. Lett. 38, 1440 (1977a).
- Peccei and Quinn (1977b) R. D. Peccei and H. R. Quinn, Phys. Rev. D16, 1791 (1977b).
- Arias et al. (2012) P. Arias, D. Cadamuro, M. Goodsell, J. Jaeckel, J. Redondo, and A. Ringwald, JCAP 1206, 013 (2012), arXiv:1201.5902 [hep-ph] .
- Svrcek and Witten (2006) P. Svrcek and E. Witten, JHEP 06, 051 (2006), arXiv:hep-th/0605206 [hep-th] .
- Arvanitaki et al. (2010) A. Arvanitaki, S. Dimopoulos, S. Dubovsky, N. Kaloper, and J. March-Russell, Phys. Rev. D81, 123530 (2010), arXiv:0905.4720 [hep-th] .
- Ringwald (2014) A. Ringwald, Proceedings, 18th International Symposium on Particles, Strings and Cosmology (PASCOS 2012): Merida, Yucatan, Mexico, June 3-8, 2012, J. Phys. Conf. Ser. 485, 012013 (2014), arXiv:1209.2299 [hep-ph] .
- Halverson et al. (2017) J. Halverson, C. Long, and P. Nath, (2017), arXiv:1703.07779 [hep-ph] .
- Kaup (1968) D. J. Kaup, Phys. Rev. 172, 1331 (1968).
- Ruffini and Bonazzola (1969) R. Ruffini and S. Bonazzola, Phys. Rev. 187, 1767 (1969).
- Das (1963) A. Das, Journal of Mathematical Physics 4, 45 (1963), http://dx.doi.org/10.1063/1.1703887 .
- Feinblum and McKinley (1968) D. A. Feinblum and W. A. McKinley, Phys. Rev. 168, 1445 (1968).
- Teixeira et al. (1975) A. F. D. F. Teixeira, I. Wolk, and M. M. Som, Phys. Rev. D12, 319 (1975).
- Colpi et al. (1986) M. Colpi, S. L. Shapiro, and I. Wasserman, Phys. Rev. Lett. 57, 2485 (1986).
- Seidel and Suen (1991) E. Seidel and W. M. Suen, Phys. Rev. Lett. 66, 1659 (1991).
- Tkachev (1991) I. I. Tkachev, Phys. Lett. B261, 289 (1991).
- Kolb and Tkachev (1993) E. W. Kolb and I. I. Tkachev, Phys. Rev. Lett. 71, 3051 (1993), arXiv:hep-ph/9303313 [hep-ph] .
- Kolb and Tkachev (1994) E. W. Kolb and I. I. Tkachev, Phys. Rev. D49, 5040 (1994), arXiv:astro-ph/9311037 [astro-ph] .
- Braaten et al. (2016a) E. Braaten, A. Mohapatra, and H. Zhang, Phys. Rev. Lett. 117, 121801 (2016a), arXiv:1512.00108 [hep-ph] .
- Eby et al. (2016a) J. Eby, P. Suranyi, and L. C. R. Wijewardhana, Mod. Phys. Lett. A31, 1650090 (2016a), arXiv:1512.01709 [hep-ph] .
- Eby et al. (2016b) J. Eby, M. Leembruggen, P. Suranyi, and L. C. R. Wijewardhana, JHEP 12, 066 (2016b), arXiv:1608.06911 [astro-ph.CO] .
- Levkov et al. (2017) D. G. Levkov, A. G. Panin, and I. I. Tkachev, Phys. Rev. Lett. 118, 011301 (2017), arXiv:1609.03611 [astro-ph.CO] .
- Helfer et al. (2016) T. Helfer, D. J. E. Marsh, K. Clough, M. Fairbairn, E. A. Lim, and R. Becerril, (2016), arXiv:1609.04724 [astro-ph.CO] .
- Braaten et al. (2017) E. Braaten, A. Mohapatra, and H. Zhang, Phys. Rev. D96, 031901 (2017), arXiv:1609.05182 [hep-ph] .
- Braaten et al. (2016b) E. Braaten, A. Mohapatra, and H. Zhang, Phys. Rev. D94, 076004 (2016b), arXiv:1604.00669 [hep-ph] .
- Bai et al. (2016) Y. Bai, V. Barger, and J. Berger, JHEP 12, 127 (2016), arXiv:1612.00438 [hep-ph] .
- Eby et al. (2017) J. Eby, M. Leembruggen, J. Leeney, P. Suranyi, and L. C. R. Wijewardhana, (2017), arXiv:1701.01476 [astro-ph.CO] .
- Desjacques et al. (2017) V. Desjacques, A. Kehagias, and A. Riotto, (2017), arXiv:1709.07946 [astro-ph.CO] .
- Dashen et al. (1975) R. F. Dashen, B. Hasslacher, and A. Neveu, Phys. Rev. D 11, 3424 (1975).
- Kudryavtsev (1975) A. E. Kudryavtsev, Soviet Journal of Experimental and Theoretical Physics Letters 22, 82 (1975).
- Makhankov (1978) V. Makhankov, Physics Reports 35, 1 (1978).
- Geicke (1984) J. Geicke, Physica Scripta 29, 431 (1984).
- Zakharov and Kuznetsov (1986) V. E. Zakharov and E. A. Kuznetsov, Zhurnal Eksperimentalnoi i Teoreticheskoi Fiziki 91, 1310 (1986).
- Gleiser (1994) M. Gleiser, Phys. Rev. D 49, 2978 (1994).
- Copeland et al. (1995) E. J. Copeland, M. Gleiser, and H. R. Muller, Phys. Rev. D52, 1920 (1995), arXiv:hep-ph/9503217 [hep-ph] .
- Salmi and Hindmarsh (2012) P. Salmi and M. Hindmarsh, Phys. Rev. D85, 085033 (2012), arXiv:1201.1934 [hep-th] .
- Mukaida et al. (2017) K. Mukaida, M. Takimoto, and M. Yamada, JHEP 03, 122 (2017), arXiv:1612.07750 [hep-ph] .
- Hormuzdiar and Hsu (1999) J. N. Hormuzdiar and S. D. H. Hsu, (1999), arXiv:hep-th/9906058 [hep-th] .
- Matos and Ureña-López (2000) T. Matos and L. A. Ureña-López, Classical and Quantum Gravity 17, L75 (2000).
- Urena-Lopez (2002) L. A. Urena-Lopez, Class. Quant. Grav. 19, 2617 (2002), arXiv:gr-qc/0104093 [gr-qc] .
- Ablowitz et al. (1973) M. J. Ablowitz, D. J. Kaup, A. C. Newell, and H. Segur, Phys. Rev. Lett. 30, 1262 (1973).
- Di Vecchia and Veneziano (1980) P. Di Vecchia and G. Veneziano, Nucl. Phys. B171, 253 (1980).
- Grilli di Cortona et al. (2016) G. Grilli di Cortona, E. Hardy, J. Pardo Vega, and G. Villadoro, JHEP 01, 034 (2016), arXiv:1511.02867 [hep-ph] .
- Petreczky et al. (2016) P. Petreczky, H.-P. Schadler, and S. Sharma, Phys. Lett. B762, 498 (2016), arXiv:1606.03145 [hep-lat] .
- Borsanyi et al. (2016) S. Borsanyi et al., Nature 539, 69 (2016), arXiv:1606.07494 [hep-lat] .
- Guth et al. (2015) A. H. Guth, M. P. Hertzberg, and C. Prescod-Weinstein, Phys. Rev. D92, 103513 (2015), arXiv:1412.5930 [astro-ph.CO] .
- Widrow and Kaiser (1993) L. M. Widrow and N. Kaiser, Astrophys. J. 416, L71 (1993).
- Marsh and Pop (2015) D. J. E. Marsh and A.-R. Pop, Mon. Not. Roy. Astron. Soc. 451, 2479 (2015), arXiv:1502.03456 [astro-ph.CO] .
- Chavanis (2011) P.-H. Chavanis, Phys. Rev. D84, 043531 (2011), arXiv:1103.2050 [astro-ph.CO] .
- Chavanis and Delfini (2011) P. H. Chavanis and L. Delfini, Phys. Rev. D84, 043532 (2011), arXiv:1103.2054 [astro-ph.CO] .
- Cotner (2016) E. Cotner, Phys. Rev. D94, 063503 (2016), arXiv:1608.00547 [astro-ph.CO] .
- Davidson and Schwetz (2016) S. Davidson and T. Schwetz, Phys. Rev. D93, 123509 (2016), arXiv:1603.04249 [astro-ph.CO] .
- Hertzberg et al. (2008) M. P. Hertzberg, M. Tegmark, and F. Wilczek, Phys. Rev. D78, 083507 (2008), arXiv:0807.1726 [astro-ph] .
- Visinelli and Gondolo (2009) L. Visinelli and P. Gondolo, Phys. Rev. D80, 035024 (2009), arXiv:0903.4377 [astro-ph.CO] .
- Derrick (1964) G. H. Derrick, Journal of Mathematical Physics 5, 1252 (1964), http://dx.doi.org/10.1063/1.1704233 .
- Rosen (1968) G. Rosen, Journal of Mathematical Physics 9, 999 (1968), http://dx.doi.org/10.1063/1.1664694 .
- Lee (1976) T. D. Lee, Phys. Rep. 23 (1976) 254-258, In *Gervais, J.L. (Ed.), Jacob, M. (Ed.): Non-linear and Collective Phenomena In Quantum Physics*, 21-25, Phys. Rept. 23, 254 (1976).
- Friedberg et al. (1976a) R. Friedberg, T. D. Lee, and A. Sirlin, Nucl. Phys. B115, 1 (1976a).
- Friedberg et al. (1976b) R. Friedberg, T. D. Lee, and A. Sirlin, Nucl. Phys. B115, 32 (1976b).
- Col (1986) Nuclear Physics B 269, 744 (1986).
- Bogolyubskiǐ and Makhan’kov (1976) I. L. Bogolyubskiǐ and V. G. Makhan’kov, ZhETF Pisma Redaktsiiu 24, 15 (1976).
- Gleiser and Sornborger (2000) M. Gleiser and A. Sornborger, Phys. Rev. E62, 1368 (2000), arXiv:patt-sol/9909002 [patt-sol] .
- Fodor et al. (2006) G. Fodor, P. Forgacs, P. Grandclement, and I. Racz, Phys. Rev. D74, 124003 (2006), arXiv:hep-th/0609023 [hep-th] .
- Piette and Zakrzewski (1998) B. Piette and W. J. Zakrzewski, Nonlinearity 11, 1103 (1998).
- Alfimov et al. (2000) G. L. Alfimov, W. A. B. Evans, and L. Vázquez, Nonlinearity 13, 1657 (2000).
- Gleiser and Sicilia (2009) M. Gleiser and D. Sicilia, Phys. Rev. D80, 125037 (2009), arXiv:0910.5922 [hep-th] .
- Gleiser and Sicilia (2008) M. Gleiser and D. Sicilia, Phys. Rev. Lett. 101, 011602 (2008), arXiv:0804.0791 [hep-th] .
- Schiappacasse and Hertzberg (2017) E. D. Schiappacasse and M. P. Hertzberg, (2017), arXiv:1710.04729 [hep-ph] .
- Chavanis (2017) P.-H. Chavanis, (2017), arXiv:1710.06268 [gr-qc] .