Dark matter bound states via emission of scalar mediators

If dark matter (DM) couples to a force carrier that is much lighter than itself, then it may form bound states in the early universe and inside haloes. While bound-state formation via vector emission is known to be efficient and have a variety of phenomenological implications, the capture via scalar emission typically requires larger couplings and is relevant to more limited parameter space, due to cancellations in the radiative amplitude. However, this result takes into account only the trilinear DM-DM-mediator coupling. Theories with scalar mediators include also a scalar potential, whose couplings may participate in the radiative transitions. We compute the contributions of these couplings to the radiative capture, and determine the parameter space in which they are important.


Introduction
In a variety of theories, dark matter (DM) is hypothesised to couple to scalar mediators. Such a coupling may determine the DM density via thermal freeze-out, and/or connect the dark sector to the Standard Model (SM) particles. It is a particularly compelling possibility in view of the discovery of the Higgs boson, which may itself be the mediator, or provide a portal to DM via an extended scalar sector. If DM possesses a sizeable coupling to the SM Higgs, then it is constrained by current experiments to be significantly heavier than the Higgs. On the other hand, if DM couples to a non-SM scalar mediator that is much lighter than itself, then this coupling may result in significant DM self-scattering inside haloes that can affect the galactic structure and bring theoretical predictions in better agreement

JHEP01(2019)070
with observations [1]. Even outside the self-interacting DM regime, the coupling of DM to a lighter non-SM scalar mediator is a generic possibility within dark sector models. In scenarios where DM couples directly to a light force carrier, non-perturbative effects associated with the long-range nature of the interaction impact the DM phenomenology. It is well known that the Sommerfeld effect [2,3] can influence the DM annihilation and self-scattering rates. This can be the case either if the mediator is a non-SM scalar or the SM Higgs [4]. More recently, it has been realised that the cosmological and astrophysical formation of DM bound states is a generic implication of theories with light mediators. The formation and subsequent decay of unstable bound states can deplete the DM density [5], and contribute to the DM indirect detection signals [6][7][8][9][10][11][12][13][14][15]. The formation of stable bound states may quell the DM self-scattering inside haloes [16], and give rise to novel radiative [17][18][19][20] and direct detection signatures [21,22]. Attractive interactions may result in the formation of large bound states [23][24][25] or non-topological solitons [26][27][28][29].
While the importance of bound-state formation (BSF) for DM is now well established for weakly coupled theories with vector mediators [5-13, 16-18, 30-39], and of course in confining DM theories [40][41][42][43][44], in models with scalar mediators bound-state effects appear to be either less severe, or relevant to a limited parameter space. In such models, the radiative capture into bound states suffers from two cancellations: (i) The lowest order s-wave contribution vanishes due to the orthogonality of the incoming and outgoing wavefunctions, thus abdicating the leading order to the p-wave. (ii) The leading order p-wave terms cancel for particle-antiparticle or identical-particle pairs, yielding to s-and d-wave contributions that are suppressed by higher orders in the coupling with respect to the capture via vector emission and possibly to the annihilation cross-section [9,24,31]. Even so, BSF has some important implications. Since it contains an s-wave component, it results in CMB constraints for fermionic DM, whose direct annihilation is p-wave and thus unconstrained by indirect probes [15]. Moreover, asymmetric DM coupled to a light scalar may form stable multiparticle bound states, provided that the coupling is large enough to overcome the two-particle capture bottleneck [24,25,45] (or BSF is facilitated by a parametric resonance [46]).
The above results take into account only the trilinear DM-DM-mediator coupling. This coupling is responsible for the long-range interaction between the DM particles that determines the scattering and bound state wavefunctions. The same coupling contributes to the radiative part of the capture process. However, there may be also other couplings that contribute to this piece. Since the radiative vertex enters the overlap integral between the initial and final states, its strength and momentum dependence are critical in determining the efficacy of BSF. It is thus important to consider all relevant contributions, particularly if they are generic within a theory. Theories with scalar mediators include a scalar potential whose couplings can contribute to the radiative part of the capture process. In this work, we investigate the impact of the various couplings in the scalar potential to the radiative BSF in theories with scalar mediators.
The paper is organised as follows. In section 2, we introduce the interaction Lagrangians, review the computation of radiative transitions, and summarise the past results on BSF with scalar emission via the trilinear DM-DM-mediator coupling. In section 3, we compute the contributions to the radiative BSF from other scalar couplings. We consider both scalar and fermionic DM, and compute BSF via one and two scalar emission (BSF1 and BSF2 respectively). We discuss the features of the resulting cross-sections and compare them with the past results. We conclude in section 4, with a discussion of their potential implications. Various technical computations are included in the appendices. For easy reference, in table 1 we summarise the notation used throughout the paper.

Preliminaries
We consider two particles X 1 and X 2 with masses m 1 and m 2 respectively, that interact via a light scalar mediator ϕ of mass m ϕ . We shall allow X 1 and X 2 to be real or complex scalars, or Dirac fermions. The relevant interactions are described by the follow-

JHEP01(2019)070
ing Lagrangians, and Note that, since we are interested in the application of our results to DM, we shall assume that the interacting particles X 1 , X 2 carry a Z 2 symmetry. For later convenience, we define the total and the reduced mass of the two interacting particles 4) and the dimensionless factors In the non-relativistic regime, the interaction between X 1 and X 2 is described to leading order by a static Yukawa potential that arises from the resummation of the oneboson-exchange diagrams, with α = α sc or α = α f , depending on whether the interacting particles are scalars or fermions, where α sc ≡ g 1 g 2 16π and We derive the Yukawa potential and α sc , α f in appendix A. The long-range interaction between X 1 and X 2 described by the potential (2.6) distorts the wavefunction of the scattering (unbound) states -a phenomenon known as the Sommerfeld effect [2, 3] -and gives rise to bound states. Bound states exist if the mediator is sufficiently light. For the ground state to exist, µα/m ϕ > 0.84 , (2.8) while stronger conditions apply for excited states [9]. The condition (2.8) also roughly marks the regime where the Sommerfeld effect is significant.

JHEP01(2019)070
The capture into bound states necessitates the dissipation of the binding energy and the kinetic energy of the relative motion of the X 1 X 2 pair, which may occur radiatively. In section 2.2, we review the computation of radiative BSF amplitudes. The capture via emission of one scalar mediator, X 1 + X 2 → B(X 1 X 2 ) + ϕ, has been previously considered in refs. [9,15,24,31], where only the trilinear ϕX 2 j couplings were taken into account. We review the main results in section 2.3. As we shall see, for a particle-antiparticle pair or a pair of identical particles, the dipole contribution -which is the leading order term for X 1 , X 2 with different masses and couplings -vanishes identically.

Radiative transition amplitude
We consider the radiative transitions where the parentheses denote the 4-momenta of the incoming and outgoing X 1 and X 2 fields. We will be interested in particular in the case where the incoming X 1 , X 2 particles form a scattering state, while the outgoing X 1 , X 2 are captured into a bound state. In order to separate the motion of the CM from the relative motion, we make the following transformation in the momenta [31,47] where η 1,2 are defined in eq. (2.5).
In the presence of a long-range interaction, the relative motion of X 1 , X 2 is not well approximated by a plane wave. It is described more generally by wavefunctions, which in momentum space we shall denote asφ k (q) andψ n m (p) for the scattering and the bound states respectively. The wavefunctions obey the Schrödinger equation with the potential (2.6). The continuous spectrum is characterised by the momentum k = µv rel , which is the expectation value of q and parametrises the energy of the relative motion in the scattering states, k = k 2 /(2µ) = µv 2 rel /2. The bound states are characterised by the standard discrete principal and angular momentum quantum numbers {n m}, which determine the expectation value of p and the binding energy n . As is well known, for a Coulomb potential, n = κ 2 /(2n 2 µ) = µα 2 /(2n 2 ), where κ ≡ µα is the Bohr momentum; a non-negligible mediator mass suppresses n and introduces a dependence on . We review the wavefunctions in appendix B (see ref. [9] for a more detailed discussion). For the purpose of evaluating the leading order contributions to the transition amplitude, we shall keep in mind that the wavefunctions impose |q| ∼ |k| = µv rel and |p| ∼ κ = µα.
In the non-relativistic regime, the total 4-momenta of the scattering and the bound states are A T Figure 1. The amplitude for the radiative capture into bound states via emission of one or two scalars consists of the initial and final state wavefunctions, and the perturbative radiative amplitude M T that includes the radiative vertices.
We will work in the CM frame, K = 0. Then, taking into account that k , n M (or equivalently α, v rel 1), the total energy available to be dissipated is (2.14) Evidently, the bound state acquires momentum |P| ∼ ω. The full amplitude for the radiative capture into a bound state depends on the overlap of the initial (scattering) and final (bound) state wavefunctions and the radiative vertex. The diagrammatic representation for transitions with emission of one or two scalars is shown in figure 1. In the instantaneous and non-relativistic approximations, the amplitude is [31] M k→{n m} 1 √ 2µ 15) where A T (q, p) is the radiative amplitude for the (off-shell) transition (2.9), under the transformations of eqs. (2.10) and (2.11). The fully connected diagrams contributing to A T can be evaluated at leading order by setting the incoming and outgoing X 1 , X 2 onshell. For any non-fully-connected diagrams contributing to A T , the virtuality of X 1 and X 2 has to be integrated out. This can be done starting from the off-shell amplitude as described in [31] (see also [34, section 2.3] for a brief summary), or by adopting an effective field theory approach [48][49][50][51][52][53][54][55][56][57][58][59]. 1 The dominant contributions to the capture with emission of a vector or scalar mediator via the trilinear coupling arise from non-fully-connected diagrams [9,31]. (However, in non-Abelian theories, a leading order contribution to capture via gluon emission arises also from a fully connected diagram [8].) In the computations of this paper in section 3, we will consider only fully connected diagrams.
In eq. (2.15), the factor inside the square brackets includes the leading order corrections in p 2 , q 2 arising from the relativistic normalisation of states [31]. Upon the convolution with the wavefunctions and integration over p and q, these terms amount to corrections in α 2 and v 2 rel . They become important when the leading order contribution from A T alone cancels, as is the case for the capture of particle-antiparticle or identical-particle pairs with emission of a scalar mediator via the trilinear coupling.

Capture with scalar emission via the trilinear coupling
The leading contributions to capture with emission of a scalar via the trilinear coupling only are shown in figure 2. Starting from eq. (2.15), and neglecting the correction arising from the normalisation of states, these diagrams yield 2 [9, 31] where P ϕ is the momentum of the emitted boson, with P 2 ϕ + m 2 ϕ = ω 2 and ω given in eq. (2.14). The wavefunctions imply that the integrand is significant in the region r 1/ max(µα/n, µv rel ) 1/ω, therefore we may evaluate eq. (2.16) by expanding in P ϕ · r ∼ max(α/n, v rel ) (cf. ref. [9, appendix B]). Clearly, the zeroth order terms vanish due to the orthogonality of the wavefunctions.
For a pair of particles with different masses and/or couplings to the scalar mediator, the leading order contributions arise from the P ϕ · r terms in the expansion of eq. (2.16). The cross-section for capture into the ground state {n m} = {100} is [9,31] where S TC BSF1 depends on the dimensionless parameters α/v rel and µα/m ϕ . In the Coulomb approximation, which is valid at m ϕ µv rel , S TC BSF1 depends only on the parameter ζ ≡ α/v rel as follows [9,31] For the more general case that includes the effect of the mediator mass, we refer to [9]. We note that eqs. (2.17) and (2.18) correspond to dipole emission ( S = 1 mode of the scattering state wavefunction). In eq. (2.17), the factor in the square brackets reduces JHEP01(2019)070 to 1 for a pair of scalar particles with g 1 = g 2 and η 1 η 2 . In eq. (2.18), the first factor inside the brackets is responsible for the characteristic σv rel ∝ 1/v rel scaling of the Sommerfeld-enhanced processes at low velocities (ζ 1), while the remaining factors tend to a constant. Evidently, for a particle-antiparticle pair or a pair of identical particles (g 1 = g 2 = g and η 1 = η 2 = 1/2), the two terms proportional to P ϕ in the expansion of eq. (2.16) cancel. This cancellation persists even for capture into bound states of non-zero angular momentum ( > 0), as is evident from eq. (2.16), and implies that the next order contributions should be considered.
For a scalar particle-antiparticle pair the next order terms, which include also the leading order correction from the normalisation of states shown in eq. (2.15), yield [9] σ TC BSF1 v rel where in the Coulomb regime S TC BSF1,XX * depends on ζ = α sc /v rel as follows (For a pair of identical scalars, an extra factor of 2 arises from the symmetrization of the scattering state wavefunction.) The case of fermionic DM was considered in ref. [24]. In the Coulomb regime and for ζ 1, the cross-section was found to be Note that both eq. (2.19) and eq. (2.21) are suppressed by α 2 with respect to eq. (2.17).
3 The contribution of the scalar potential couplings to radiative capture In this section, we investigate how the scalar couplings λ 1ϕ , λ 2ϕ , λ ϕ and ρ ϕ in eqs. (2.1) to (2.3) -which do not affect the long-range potential between X 1 and X 2 -contribute to the radiative capture of X 1 , X 2 pairs into bound states. In order to exhibit the leading order contribution from all couplings, we consider BSF via one and two scalar emission, It is important to note that in most diagrams we will consider, the trilinear couplings g j also participate in the radiative part of the process (cf. figures 3 to 6). In fact, it may naively seem that these diagrams are of the same or higher order in α than those giving rise to the cross-sections of section 2.3 (even ignoring the additional suppression introduced by the new couplings). However, the momentum transfer along the mediators exchanged in these diagrams and the off-shellness of the interacting particles scale also with α, thereby reducing the order of dependence of the diagrams on α. 3

JHEP01(2019)070
For convenience, we first define in section 3.1 the wavefunction overlaps integrals that we will use for the BSF cross-sections. We provide some analytic approximations for capture into the ground state, and describe their features that are of course inherited by the BSF cross-sections. Then, in sections 3.2 and 3.3 we compute the contributions to BSF1 from the couplings in the scalar potential, for scalar and fermionic interacting particles, respectively, and discuss in which regimes they are important. We finish with considering BSF2 in section 3.4.

Overlap integrals
We define the overlap integrals The prefactors in eqs. ( where we used the Hulthen approximation for the binding energy (cf. appendix B). In the regime where BSF is important, v rel α (see e.g. section 2.3), this condition reduces roughly to m ϕ µα 2 /2 or equivalently ξ ξ min 2.4/α 1, i.e. it is much stronger than the requirement for the existence of bound states. This in turn ensures that the bound-state wavefunction can be approximated by its Coulomb value. (For the validity of the Coulomb limit for the overlap integrals and the BSF cross-sections, see below. ) We derive analytical expressions for the V k,{100} and R k,{100} integrals in appendix C, using an appropriate approximation. The integral (3.5) has been considered in ref. [9]. For the BSF cross-sections of interest, |Γ| ∼ |P ϕ | |p| ∼ κ, thus I k→{n m} can be computed by expanding in Γ. Here we shall keep up to first order terms in Γ/κ. For capture into the JHEP01(2019)070 ground state, where S 0 (ζ, ξ) and S 1 (ζ, ξ) are the Sommerfeld factors for s-and p-wave annihilation respectively. In the Hulthen approximation, . While there is no analytical approximation for S 1 (ζ, ξ) for finite ξ, in the Coulomb limit ξ → ∞ it becomes S C 1 (ζ) = (1 + ζ 2 )S C 0 (ζ). The Coulomb limit remains a good approximation as long the average momentum transfer between the interacting particles is larger than the mediator mass, m ϕ µv rel , (3.12) or equivalently ξ > ζ (see e.g. [9]). Outside this range, i.e. at low velocities v rel m ϕ /µ, both S 0 and S 1 exhibit parametric resonances at discrete values of ξ that correspond to the thresholds for the existence of = 0 and = 1 bound states respectively. For S 0 in the Hulthen approximation, these are ξ = n 2 with n ∈ integers. At non-resonant parametric points, S 0 and S 1 follow the Coulomb approximation as long as (3.12) is satisfied, but saturate to their respective Coulomb values at v rel ≈ m ϕ /µ as the velocity decreases. In contrast, if close to a resonant parametric point, S 0 and S 1 grow faster than 1/v rel and 1/v 3 rel respectively at v rel m ϕ /µ (in particular S 0 ∝ ζ 2 ∝ 1/v 2 rel ), to eventually saturate to a constant value at a lower velocity that depends on the proximity to the resonant point. For S 0 in the Hulthen approximation, the saturated value is S 0 (ζ, ξ) π 2 ξ/ sin 2 π 2 ξ [cf. eq. (3.11)].
This behaviour implies that an s-wave cross-section times relative velocity saturates to a constant value at low velocities, while a p-wave recovers the velocity suppression that appears in perturbative cross-sections, σ p−wave v rel ∝ v 2 rel , albeit is enhanced with respect to its value if the Sommerfeld effect were neglected. This point will be important in our discussion of the new contributions to BSF that we compute in the following.
3.2 Capture via emission of one scalar mediator ϕ, scalar X 1 and X 2

Amplitude
The leading order diagrams in λ 1ϕ , λ 2ϕ and ρ ϕ that contribute to the radiative part of the transition amplitude due to the scalar couplings in eqs.  Figure 3. The contribution of the scalar couplings to the radiative part of bound-state formation via emission of one scalar mediator. For fermionic X 1 , X 2 , only the diagram to the right exists.

JHEP01(2019)070
We find (3.13) Then according to the discussion in section 2.2 on the scaling of the momenta, (3.14) The contribution from the diagrams of figure 2 is [9,31] A TC

Cross-section for capture into the ground state
The cross-section for capture into the ground state is The momentum of the emitted scalar is found from the conservation of energy [cf. eq. (2.14)] where we adopted the Coulomb value for the binding energy since ξ 1.

JHEP01(2019)070
Taking into account the contributions of eqs. (3.16) and (3.17) to the amplitude, and the overlap integrals of eqs. (3.8) to (3.10) for capture into the ground state, we obtain where s BSF1 ps is the phase-space suppression factor for BSF1, As can be seen from eq. (3.20), the ρ ϕ contribution to the BSF1 cross-section scales as σ BSF1 v rel ∝ (1/α)(ρ ϕ /µ) 2 . This is clearly of lower order in α than the cross-sections of section 2.3, however it is suppressed by the square of the ratio of the mediator scale to the DM scale. Assuming that ρ ϕ ∼ m ϕ , the kinematic threshold m ϕ µα 2 /2 implies that this contribution scales at best as σ BSF1 v rel ∝ α 3 . While this is of higher order in α than eq. (2.17), it is still of lower order than eqs. (2.19) and (2.21). It is of course important to keep in mind that ρ ϕ may differ significantly from m ϕ and/or the binding energy. 5 We shall now simplify and adapt the above expression -which has been derived for a pair of distinguishable scalars -to various cases.
Different scalars species. Let us consider for simplicity the limiting case where the couplings of the two particles are equal, λ 1ϕ = λ 2ϕ = λ Xϕ and g 1 = g 2 = g, while their masses are very different, η 1 η 2 . Then, eq. (3.20) simplifies to . (3.22) In the regime where BSF is important and the Coulomb approximation is valid, i.e. for m ϕ µv rel µα sc , we recall that S 1 (1 + ζ 2 )S 0 , thus all contributions exhibit the same velocity scaling, σv rel ∝ 1/v rel . For perturbative λ Xϕ and for ρ ϕ ∼ m ϕ µα 2 sc /2, the s-wave contributions are subdominant, and we recover the cross-section arising from the diagrams of figure 2 [cf. eq. (2.17)]. 5 We note that the third diagram of figure 3, from where the ρϕ contribution to BSF1 arises, resembles the diagram that appears in the radiative capture via one gluon emission in non-Abelian theories due to the trilinear gluon coupling [8]. That coupling is momentum dependent, with the relevant momentum scale in the capture process being κ. We may recover the scaling of the non-Abelian BSF cross-section on α [34] by mapping ρϕ → g κ ∝ α 3/2 .

JHEP01(2019)070
However, outside the Coulomb regime, i.e. for µv rel m ϕ µα sc , the p-wave term becomes velocity suppressed (cf. discussion in the end of section 3.1). This implies that at sufficiently low velocities, the λ Xϕ and ρ ϕ contributions dominate.
Note that if the mediator ϕ is the radial component of a complex scalar Φ that obtains a vacuum expectation value v ϕ and breaks a local symmetry, i.e. Φ = (v ϕ +ϕ)e iaϕ / √ 2, then the trilinear DM-DM-mediator coupling g arises from the quartic coupling of this scalar to DM after spontaneous symmetry breaking, i.e.
In this case, λ Xϕ and g are related via λ Xϕ v ϕ = gm X . Note that m X receives a contribution from v ϕ , but remains an independent parameter. The above implies α sc = λ 2 Xϕ v 2 ϕ /(16πm 2 X ), and the condition (3.24) becomes λ Xϕ (16π/18)(m X /v ϕ ) 2 , which encompasses the entire regime where λ Xϕ is perturbative if v ϕ ∼ m ϕ m X (or more generally, if v ϕ m X ). In this case, the λ Xϕ contribution to BSF1 dominates.

JHEP01(2019)070
Identical scalars. For a pair of identical scalars, the total wavefunction has to be symmetric in the interchange of the two particles. Thus the scattering state wavefunction is

27)
This implies that the contribution of the even-S modes participating in a process is doubled with respect to the case of distinguishable scalars, while the contribution of the odd-S modes vanishes. Similarly, there are only = even bound states of two identical bosons. Since eq. (3.23) includes only s-wave terms, the cross-section for a pair of identical scalars is twice as large as that given by eq. (3.23).
3.3 Capture via emission of one scalar mediator ϕ, fermionic X 1 and X 2

Amplitude
For fermions, there are no renormalisable λ jϕ couplings and only the third diagram in figure 3 contributes. The amplitude is related to that for scalars as follows. The Dirac propagators S D j can be expressed in terms of the scalar propagators S j (p) and the spinors u,ū where r denotes the spin and j = 1, 2 refers to the particle species. In order to compute the fermionic BSF diagrams, we insert a factor r i u r i j (p i )ū r i j (p i ) for each propagator and contractū r i+1 j (p i+1 )u r i j (p i ) across each vertex i. Since all the fermion-fermion-scalar vertices in the BSF diagrams are either soft or ultrasoft, we can use the identitȳ as also in appendix A. The spin Kronecker deltas in eq. (3.29), upon summation over the internal spin indices, ensure that the spin of each particle is conserved across the entire diagram, including both the ladders and the vertices in the radiative parts of the diagrams. With this, we find where r 1 , r 1 and r 2 , r 2 are the spins of the incoming and outgoing X 1 and X 2 particles. Comparing eq. (3.30) to (3.14), we see that the ρ ϕ contribution to A T is larger by a factor 4 for a fermionic pair than for a scalar pair. However, taking into account that for fermions α f = g 1 g 2 /(4π), the full amplitude for fermions looks the same as the ρ ϕ contribution to eq. (3.16), up to the spin conservation factors,

JHEP01(2019)070
where we used eqs. (2.15) and (3.4). Similarly, the contribution to the amplitude from the diagrams of figure 2 is where the extra factor 2 with respect to eq. (3.17) arises from the spinor contraction along the leg that contains the radiative vertex.

Cross-section for capture into the ground state
Upon squaring the amplitudes (3.31) and (3.32), summing over the final-state spins and averaging over the spin of the initial particles, the spin factors simply yield 1. Using eqs. (3.18) and (3.19) and the overlap integrals (3.9) and (3.10), we find the spin-averaged BSF1 crosssection to be As in the previous section, we shall now consider some specific cases.
Particle-antiparticle pair. For a particle-antiparticle pair, the p-wave contribution in eq. (3.33) vanishes, and the BSF cross-section becomes We recall that there are also s-and d-wave contributions of higher order in α f from the trilinear ϕXX coupling alone. Comparing eq. (3.35) with eq. (2.21), we find that the ρ ϕ term is more significant if

JHEP01(2019)070
Identical fermions. For a pair of identical fermions, the total wavefunction has to be antisymmetric in the interchange of the two particles. This implies that the spatial wavefunction depends on their total spin. A pair of identical spin-1/2 particles may be either in the antisymmetric spin-0 state, or in the symmetric spin-1 state. Their spatial wavefunction should then be symmetric or antisymmetric, respectively, fermions with total spin 0: fermions with total spin 1: As for identical bosons, the wavefunction (3.37) implies that the contribution of the even-S modes participating in a process is doubled with respect to the case of distinguishable particles, while the contribution of the odd-S modes vanishes. The opposite holds for the wavefunction (3.38). Similarly, there are only = even bound states of two identical totalspin-0 fermions, and only = odd bound states of two identical total-spin-one fermions. Since eq. (3.35) is an s-wave process, the spin-averaged cross-section for capture of two identical fermions is half of that given in eq. (3.35), with the entire contribution arising from the spin-0 state.

Capture via emission of two scalar mediators
The scalar couplings appearing in eqs. (2.1) to (2.3) raise the possibility that the radiative part of the capture into bound states may carry a lower suppression in α if two mediators are emitted. On the other hand, the emission of two mediators, which share the available energy ω µ(α 2 + v 2 rel )/2, implies that the BSF2 cross-section picks up additional suppression in powers of α (for α > v rel ) due the phase space of the second mediator. To determine which diagrams may contribute significantly, we first work out the phase-space integration. 6

Phase-space integration
The cross-section for capture via emission of two mediators with momenta P a and P b is (3.39) where the factor (1/2) is due to the two identical particles in the final state. We work in the CM frame, K = 0, and use the three-momentum delta function to integrate over the momentum of the bound state P = −(P a + P b ). The energy delta function yields the condition In a previous version of this paper, it was stated that the contributions to BSF2 from the trilinear DM-DM-mediator couplings gj alone do not suffer from the cancellations of the gj contributions to BSF1 shown in figure 2. However, the resummation of the one-scalar-exchange diagrams between the two emission vertices in the corresponding BSF2 diagrams reveals that the same cancellations are present. Thus, these contributions to BSF2 are not of interest, and we do not consider them here. We thank Joan Soto for pointing this out.

JHEP01(2019)070
For convenience, we define the dimensioneless parameters Then, combining energy and momentum conservation, we obtain the phase space condition for the momenta of the emitted scalars, where θ ab is the angle between P a and P b . Because δ 1, the phase-space encompassed by eq. (3.43) extends on the x a − x b plane essentially along the line with 0 x a √ 1 − 2d, and a small width along the x b direction, due to cos θ ab ranging in [−1, 1]. Note that for the capture process to be kinematically allowed, 2d < 1. The x b width can be estimated by differentiating eq. (3.43) with respect to x b and cos θ ab . We find For the diagrams of interest, the amplitude is independent of P a,b (cf. section 3.4.2). Thus, putting everything together, eq. (3.39) yields where s BSF2 ps is the phase-space suppression factor due to the kinematic threshold for BSF2, with s BSF2 ps = 1 for d = 0. Comparing eq. (3.46) with eq. (3.18), we observe that the BSF2 cross-section is proportional to two extra powers of the available energy to be dissipated, ω n + k , with respect to BSF1, as expected from the phase-space element for the second mediator This implies an extra suppression by α 4 in the regime where BSF is important, α > v rel . We also note the suppression of BSF2 with respect to BSF1 by the numerical factor 1/(2 4 3π 2 ) 2 × 10 −3 .

Amplitude and cross-section
In figures 4 to 6, we show various contributions to BSF2 arising from the mediator selfcouplings and the mediator-DM quartic coupling for scalar DM. The diagram of figure 4 yields the most important contribution. We shall first compute the cross-section arising from this diagram and then discuss why we neglect the other contributions. Figure 4. Capture via emission of two scalar mediators: the contribution of the mediator quartic self-coupling to the radiative amplitude.

JHEP01(2019)070
The amplitude for the λ ϕ contribution can be obtained from the second term in eq. (3.16) by replacing ρ ϕ with λ ϕ . Taking into account the discussion in section 3.3.2, we find it to be the same for both bosons and fermions, where R k,{n m} is defined in eq. (3.4). From eqs. (3.9) and (3.46) we obtain the (spinaveraged) cross-section for capture into the ground state, where s BSF2 ps and S 0 are given in eqs. (3.47) and (3.11) respectively. The factor f follows from the discussion in sections 3.2 and 3.3, with Comparing with the BSF1 contributions from the trilinear coupling alone (2.19) and (2.21), the λ ϕ contribution to BSF2 is more significant for The perturbativity of the couplings implies that the above condition may be meaningfully satisfied only for small coupling α 10 −3 .

Other subdominant contributions
In figures 5 and 6, we display various diagrams arising from ρ ϕ and λ jϕ . Their contributions to BSF2 are subdominant compared to the contributions of the same couplings to BSF1

JHEP01(2019)070
for reasonable choices of the parameters, as we discuss below. Since for phenomenological purposes, we are interested mostly in the total capture rate, we do not compute these contributions in detail. Because part of the suppression arises from the phase space density of the second mediator (cf. section 3.4.1), we will consider directly the contributions of these diagrams to the cross-section rather than the amplitude, thus neglecting any cross-terms between different diagrams (except in the cases where there is a cancellation between them). It is straightforward to verify that the cross-terms do not alter our conclusions. We focus on particle-antiparticle or identical-particle pairs. Our comparisons refer to the regime α > v rel where BSF is important, setting aside the ζ-dependent factors that yield the same v rel dependence between BSF1 and BSF2 in this regime.
• For particle-antiparticle or identical-particle pairs, the sum of the two diagrams on the left in figure 5 suffers from the same cancellation as the diagrams of figure 2 (cf. section 2.3). This introduces a α 4 suppression in the corresponding crosssection, while the phase-space density of the second emitted mediator implies an additional suppression by (µα 2 /2) 2 as discussed in section 3.4.1. These suppressions are balanced by the propagator of the off-shell mediator, which introduces a factor ∼ (P 0 a + P 0 b ) −4 ∝ (µα 2 /2) −4 . Thus, the contribution of these diagrams to the crosssection scales as σ BSF2 v rel ∝ (ρ ϕ /µ) 2 . This is suppressed by α with respect to the ρ ϕ contribution to BSF1, σ BSF1 v rel ∝ (1/α)(ρ ϕ /µ) 2 [see e.g. eq. (3.35)], as well as by the numerical factor ∼ 2 × 10 −3 due to the three-body phase space. For the third diagram, the off-shell mediator yields a factor ∝ (µα 2 /2) −4 . Taking into account the phase-space suppression ∝ (µα 2 /2) 2 , we find that the contribution to the BSF2 cross-section scales as σ BSF2 v rel ∝ (1/α 5 )(ρ ϕ /µ) 4 . Assuming that ρ ϕ ∼ m ϕ µα 2 /4 (cf. footnote 5), effectively this yields at best σ BSF2 v rel ∝ α 3 , similarly to the BSF1 cross-section (3.35). Still, the numerical suppression of BSF2 due to the three-body phase space ensures that the BSF2 contribution is subdominant to BSF1. However, if there is a large hierarchy between m ϕ and ρ ϕ with ρ ϕ m ϕ , then the diagram under consideration may be significant. 7 For the forth diagram, the off-shell mediator between the two emitted scalars yields a factor ∝ (µα) −4 . It follows that it is subdominant to the third diagram in figure 5 and to the ρ ϕ contribution to BSF1.
• For the two diagrams on the left in figure 6, we start from the first term of eq. (3.23), replace one factor of α with λ 2 Xϕ and introduce the phase-space suppression ∝ α 4 . We obtain σ BSF2 v rel ∝ α 5 λ 4 Xϕ . 7 We note that the three left diagrams in figure 5 exhibit a collinear divergence at mϕ → 0 that would have to be treated appropriately. Figure 5. Capture via emission of two scalar mediators: the contributions of the mediator trilinear self-coupling to the radiative amplitude. Figure 6. Capture via emission of two scalar mediators: the contribution of the quartic DMmediator coupling to the radiative amplitude. These diagrams exist for scalar DM only.
• The two diagrams on the right in figure 6 suffer from the cancellations of the BSF1 diagrams of figure 2. Replacing one factor of α with λ 2 Xϕ results in the suppression λ 2 Xϕ α 3 for particle-antiparticle or identical-particle pairs. Introducing also the phasespace suppression yields σ BSF2 v rel ∝ λ 2 Xϕ α 7 .

Conclusion
Dark matter coupled to a light scalar mediator is motivated in the context of Higgs portal models, self-interacting DM, and more generally models beyond the SM with extended scalar sectors. Light force mediators generically imply the existence of DM bound states, whose formation has important implications for the DM phenomenology. While the cosmological and astrophysical formation of bound states is significant in models with vector mediators, it is less efficient in models with scalars mediators. Indeed, considering only a trilinear DM-DM-mediator coupling, the radiative capture of particle-antiparticle or identical-particle pairs into bound states with emission of a scalar mediator is subject to cancellations that suppress the cross-section by higher powers in this coupling.

JHEP01(2019)070
While the trilinear coupling determines the long-range interaction between the DM particles, it is not necessarily the only contributor to the radiative part of the transition. Quite generically, in theories with scalar mediators the couplings in the scalar potential also contribute. In this work, we investigated the contribution of the mediator self-couplings, as well as the DM-mediator quartic coupling in the case of scalar DM, to the radiative capture into bound states. We considered capture both via one and two scalar mediator emission. Our main results include eqs. (3.22) and (3.23) for BSF1 by scalar particles, eqs. (3.34) and (3.35) for BSF1 by fermionic particles, and in eq. (3.49) for BSF2. We have found that the newly considered couplings can enhance or dominate the capture rate in sizeable parts of the parameter space, described in part by eqs. (3.24), (3.26), (3.36), (3.51) and (3.52).
Importantly, the new contributions are s-wave and remain significant even at very low velocities, thereby enhancing the radiative signals that can be probed in the indirect searches. This can potentially strengthen the resulting constraints. While models with light mediators and s-wave annihilation are severely constrained by the CMB and other indirect probes [10,61,62], models with p-wave annihilation -such as fermionic DM coupled to a scalar mediator -remain largely unconstrained. In this case, the formation and subsequent decay of unstable bound states offers a source of detectable indirect signals [15]. (Bremsstrahlung of dark mediators has also been invoked to lift the p-wave suppression [63].) Moreover, s-wave annihilation processes may lead to a period of reannihilation in the early universe, after DM kinetic decoupling [64].
Besides the signals produced from the decay of unstable bound states, the radiation emitted during the capture process is another source of indirect signals even in the case of stable bound states [19]. This is particularly important for asymmetric DM [65], whose annihilation signals are suppressed due to the asymmetry (albeit can still be significant due to the Sommerfeld effect [12]). Asymmetric DM, in turn, offers an excellent host of self-interacting DM, both because it allows for large couplings to light mediators, and it evades the indirect detection constraints, provided that the asymmetry is sufficiently large [13]. The enhanced radiative signals expected due to the contributions computed here can improve the prospects of probing asymmetric and self-interacting DM.
Moreover, for asymmetric DM coupled to a light scalar, the enhanced rate of formation of two-particle bound states can result in more efficient cosmological formation of multiparticle bound states [24], which in the case of scalar DM may lead to solitosynthesis [66][67][68].

A Non-relativistic potential for scalars and fermions
The non-relativistic static potential describing the X 1 − X 2 interaction is (see e.g. [31,69]) where M 2PI (k) encompasses all 2-particle-irreducible diagrams that contribute to the X 1 − X 2 elastic scattering. To leading order, this is the one boson exchange shown in figure 7. Figure 7. One boson exchange diagrams that yield the leading order contribution to the nonrelativistic potential between two different paticles X 1 and X 2 (left), or a particle-antiparticle pair (right). The parentheses denote the spin of the incoming and outgoing particles, in the case of fermions.

JHEP01(2019)070
From the Lagrangians of eqs. (2.1) to (2.3), we find for the interaction of a scalar pair X 1 X 2 , a fermionic pair X 1 X 2 and a fermion-antifermion pair XX, From eqs. (A.1) to (A.4), we find the well known Yukawa potential with α = α sc or α = α f depending on whether the interacting particles are scalars or fermions respectively, where The interaction is attractive if g 1 g 2 > 0.

JHEP01(2019)070 B Scattering state and bound state wavefunctions
The bound state and scattering state wavefunctions obey the Schrödinger equations Coulomb limit, m ϕ → 0. The ground-state and scattering state wavefunctions are where 1 F 1 is the confluent hypergeometric function of the first kind, and (B.5) φ C k (r) can also be decomposed in angular momentum modes as follows Outside the Coulomb regime: Hulthen approximation. For the Yukawa potential (2.6), the solutions of the Schrödinger equation are not known in a closed form. However, it is possible to obtain analytical solutions for the = 0 modes of the wavefunctions, if we replace the Yukawa for the Hulthen potential [70,71], For m * ∼ m ϕ , the Hulthen potential reproduces the behavior of the Yukawa potential at short and large distances. Both potentials admit bound state solutions provided that the screening length scale is sufficiently large. The thresholds for the n-th bound level in the Hulthen potential are m * 2µα/n 2 , while for the Yukawa potential (2.6), bound states exist if m ϕ µα/0.84 [9]. We shall pick m * = 1.68m ϕ , (B.8) such that the conditions for the existence of the lowest bound state coincide. The solutions to the Schrödinger equation for the Hulthen potential can be expressed in terms of the two dimensionless parameters ζ ≡ α/v rel and ξ ≡ 2µα/m * . (B.9)

JHEP01(2019)070
Note that ζ and ξ are always defined using the appropriate fermionic or scalar coupling α. ζ compares the average momentum transfer between two unbound particles (∼ µv rel ) with the relative average momentum of the particles inside the bound state (∼ µα), while ξ compares the Bohr momentum (κ = µα) that determines the size of the bound state, with the screening scale (m * ) that determines the range of the interaction. The interaction manifests as long-range roughly if ξ 1; this is the regime where non-perturbative phenomena, the Sommerfeld effect and bound states, emerge. The Coulomb limit is recovered for ξ → ∞.
For the ground state, the binding energy is which reduces to eq. (B.5) in the limit ξ → ∞.

C Overlap integrals
For the capture into bound states via emission of two scalars, we need to compute the overlap integrals In the second step in eq. (C.6), we expanded the decaying exponential inside the integral. The bound-state wavefunction implies that the integrand is significant for r n/(µα) while the kinematic threshold (3.7) for BSF imposes m ϕ < µα 2 /(2n 2 ); therefore m ϕ r α 1. The zeroth order term in the expansion vanishes due to the orthogonality of the wavefunctions, leaving the first order term to be the dominant contribution. In the following, we will focus on the capture into the ground state, {n m} = {100}.
While in the Coulomb limit (m ϕ → 0) it is possible to obtain analytical expressions for the integrals eqs. (C.5) and (C.6), outside the Coulomb regime they should be evaluated numerically. However, it is possible to obtain an analytical approximation, using the Hulthen potential (cf. appendix B). Since ξ ≡ 2µα/m * 1 in all the parameter space where BSF is kinematically allowed to occur [cf. eq. (3.7)], the bound state wavefunction ψ 100 (r) can be well approximated by its Coulomb value. On the other hand, the scattering-state wavefunction φ k (r) is close to its Coulomb limit for r 1/m * , up to an overall normalisation which, for the S = 0 mode, is determined by the factor S 0 (ζ, ξ) of eq. (B.17). Since ξ 1, this encompasses all the range in which ψ 100 (r) and therefore the integrands in eqs. (C.5) and (C.6) are important, r 1/(µα). It follows that V k,{100} and R k,{100} -which depend only on the S = 0 mode of the scattering state wavefunction -can be well approximated by their Coulomb values even outside the Coulomb regime provided that we replace S C 0 (ζ) in eq. (B.4) with S 0 (ζ, ξ).