Multicomponent Dark Matter from Gauge Symmetry

The composition of Dark Matter (DM) remains an important open question. The current data do not distinguish between single- and multi-component DM, while in theory constructions it is often assumed that DM is composed of a single field. In this work, we study a hidden sector which naturally entails multicomponent DM consisting of spin-1 and spin-0 states. This UV complete set-up is based on SU(3) hidden gauge symmetry with the minimal scalar field content to break it spontaneously. The presence of multiple DM components is a result of a residual Z_2 x Z'_2 symmetry which is part of an unbroken global U(1) x Z'_2 inherent in the Yang-Mills systems. We find that the model exhibits various parametric regimes with drastically different DM detection prospects. In particular, we find that the direct detection cross section is much suppressed in large regions of parameter space as long as the Standard Model Higgs mixes predominantly with a single scalar from the hidden sector. The resulting scattering rate is often beyond the level of sensitivity of XENON1T, while still being consistent with the thermal WIMP paradigm.


Abstract
The composition of Dark Matter (DM) remains an important open question. The current data do not distinguish between single-and multi-component DM, while in theory constructions it is often assumed that DM is composed of a single field. In this work, we study a hidden sector which naturally entails multicomponent DM consisting of spin-1 and spin-0 states. This UV complete set-up is based on SU(3) hidden gauge symmetry with the minimal scalar field content to break it spontaneously. The presence of multiple DM components is a result of a residual Z 2 × Z 2 symmetry which is part of an unbroken global U (1) × Z 2 inherent in the Yang-Mills systems. We find that the model exhibits various parametric regimes with drastically different DM detection prospects. In particular, we find that the direct detection cross section is much suppressed in large regions of parameter space as long as the Standard Model Higgs mixes predominantly with a single scalar from the hidden sector. The resulting scattering rate is often beyond the level of sensitivity of XENON1T, while still being consistent with the thermal WIMP paradigm.

Introduction
The existence of a Dark Matter (DM) component of the Universe is confirmed by several astrophysical and cosmological probes, e.g. the CMB [1] and structure formation. Particle physics solutions to the DM problem mostly rely on the existence of a new particle that is stable on cosmological scales thanks to a symmetry, and with weak enough interactions with Standard Model (SM) states to evade constraints from direct and indirect searches. From the model-building point of view, these requirements are naturally satisfied by a "hidden" sector whose states are singlets with respect to the SM symmetry group. In this kind of a setup, stable particles of the hidden sector are Dark Matter candidates.
Despite the different symmetry groups acting on the visible and hidden sectors, renormalizable interactions can arise among them. The dimension-4 operators relevant to our study are obtained by combining the gauge invariant dimension two terms H † H and B µν with similar dimension-two operators formed by the states of the hidden sector. The strength of such "portal" interactions can be sufficient to bring the visible and hidden sectors in thermal equilibrium in the Early Universe and realize the WIMP paradigm. In addition, the DM candidates retain interactions with the visible particles at present times, which could be within the reach of the current and future searches for new particles.
In the simplest models of this type, the hidden sector is populated (effectively) by a single field which constitutes DM. The lowest order Higgs portal operators mediating interactions between the DM and the SM states then read H † H|χ| 2 or H † HV µ V µ , with χ and V µ being a scalar and a vector DM candidate, respectively, see e.g. [2][3][4][5][6][7][8][9]. 1 These simple set-ups are currently under pressure from constantly improving experimental constraints. Indeed, by crossing symmetry arguments, there is a relation between the DM pair annihilation cross-section at freeze-out, responsible for the relic density, and the processes potentially responsible for detection signals, such as scattering on nuclei, probed by Direct Detection experiments, or production at colliders. Null results of the latter then rule out large portions of the parameter space favoured by the thermal WIMP paradigm [15,16]. If the s-wave DM annihilation cross section remains substantial at present times, the model can also be probed by Indirect Detection experiments. These considerations motivate exploration of richer hidden sector structures. For example, additional annihilation channels into dark sector states can deplete the DM relic density without changing its interactions with the visible sector while satisfying the experimental constraints [17,18].
More generally, there is no a priori reason for the DM of the Universe to be composed of a single field. Multi-component DM frameworks, with two or more particles contributing a non-negligible fraction to the total relic density Ω DM,tot h 2 ≈ 0.12, offer interesting perspectives. The relation between the annihilation cross section and the current detection signals has to be properly reconsidered. Some work in this direction has been carried out in [19,20], where the discovery potential of the current and future experimental facilities and the capability of discriminating multicomponent DM from single-component DM have been studied. We note that multicomponent DM emerges in various particle physics models (see e.g. [21][22][23][24]).
In this work, we will investigate multicomponent DM emerging from a hidden sector endowed with gauge symmetry. Such systems enjoy natural discrete symmetries which can act as DM stabilizers. Indeed, it was noted in [6] and detailed in [8] that a hidden sector consisting of a U(1) gauge field A µ and a single complex scalar which breaks the symmetry spontaneously has the symmetry As a result, the massive vector field A µ is stable and can constitute DM. This idea generalizes to non-Abelian gauge symmetries as well. In particular, hidden SU(N) sectors with a minimal matter content necessary to break the gauge symmetry completely, that is N − 1 scalar N -plets, are endowed with a Z 2 × Z 2 symmetry [25] 2 where A a µ is the gauge field, a is the adjoint group index and n a , n a take on values 0,1 depending on a. One of the Z 2 's corresponds to complex conjugation of the SU(N) group elements, while the other is a gauge transformation. In the SU(2) case, the symmetry enlarges in fact to custodial SO(3) [6] (see also [26,27]), while for larger groups it is part of a global unbroken U (1) × Z 2 . Since the SM fields are neutral under the above symmetries, A a µ cannot decay into the visible sector particles and therefore can constitute dark matter. Note that since the symmetry is Z 2 × Z 2 , one expects at least two different DM components, in contrast to traditional Z 2 -invariant WIMP models. The exact composition depends on the details of the spectrum. In particular, the SU(3) example studied in [25] has only vector DM with two components being degenerate in mass and the third one being somewhat lighter. These states interact with the visible sector through the Higgs portal operators. A related study has recently appeared in [28].
In our current study, we explore a qualitatively different case of mixed spin DM, that is containing both spin 1 and spin 0 components. We employ the model of [25] in a different parametric regime, where a stable pseudoscalar is lighter than the gauge field with the same Z 2 ×Z 2 quantum numbers. In this case, the pseudoscalar as well as the gauge fields with distinct Z 2 × Z 2 quantum numbers constitute DM. The resulting phenomenology is very different from that of [25]. In particular, we find that there are substantial regions of parameter space where the direct detection cross section is suppressed.
We stress that although we study a specific model of multicomponent DM, many of the results presented here are of general relevance. In particular, depending on the composition of DM, the direct detection signal strength varies drastically, over orders of magnitude, and is often consistent with thermal relic DM abundance. Such behaviour is specific to more complicated hidden sectors within our framework and reflects the possibility that common models may oversimplify the DM properties.
One of the novel aspects of our study is that multicomponent DM is a natural consequence of our UV-complete framework, due to Z 2 × Z 2 being part of the Yang-Mills symmetries. This is in contrast to more conventional models where the two DM components have different origins such as the mixed axion-neutralino DM scenario [29,30]. Consequently, the contributions of the components to the total DM density are controlled by a set of the UV parameters. In our study, much emphasis will be given to the analysis of the DM production processes (as opposed to the approach of [19,20]). We solve numerically the coupled Boltzmann equations and calculate the individual relic abundances as a function of the parameters of the model. The composition of DM can be very different in different parameter regions and in some of them both DM components give comparable contributions. We then study the Direct Detection constraints and observe an interesting effect. As long as the SM Higgs mixes predominantly with one of the hidden scalar fields, the direct detection is highly suppressed in the parameter regions where DM is mostly spin-0. It can be so small that even future detectors like XENON1T [31] will not be able to probe it. This is one of the main results of our study.
The paper is structured as follows. The model is introduced in Section 2. Section 3 is devoted to a detailed discussion of the relic DM density calculations and the Direct Detection limits. We also comment on the possibility of detecting one of the components through Indirect Detection. Our results are summarized in Section 4.

The SU(3) hidden sector model
The purpose of this section is to briefly summarise our model, mostly following Ref. [25]. The hidden sector of the model is endowed with SU(3) gauge symmetry, which is broken spontaneously (to nothing) by two hidden triplets φ 1 and φ 2 . This is the minimal setup that allows one to make all the SU(3) gauge fields massive.
The Lagrangian of the model is where Here, H is the Higgs doublet, which in the unitary gauge can be written as H T = (0, v + h)/ √ 2, and the most general renormalisable hidden sector scalar potential is given by The fields φ 1 and φ 2 are responsible for the spontaneous breaking of the hidden SU(3) symmetry. In the unitary gauge, they can be written as where the v i are real VEVs and ϕ i are real scalar fields. Here we assume the CP symmetry in the scalar sector, i.e. that the couplings are real and ϕ 4 attains no VEV. As a consequence there is no mixing between the CP-even scalar fields ϕ 1−3 and the CP-odd scalar ϕ 4 . This allows for the possibility that ϕ 4 is stable.
A minor technical complication that occurs in this model is that the quadratic part of the Lagrangian is not diagonal, due to the mixing terms These terms of the form κ ai A a µ ∂ µ ϕ i can be removed by the transformation where M is the mass matrix of the hidden gauge bosons. This leaves M unchanged. After the above transformation the kinetic terms of the ϕ i are not canonically normalised anymore so that a further transformation is needed to make the quadratic part of the Lagrangian canonically normalised. To stress its special role as a DM candidate, we relabel To simplify the analysis, in the rest of the paper we assume that the couplings m 2 12 , λ H12 , λ 6 , λ 7 in the scalar potential, as well as the VEV v 3 , are small but nonvanishing. If they did vanish, the system would attain an additional unwanted Z 2 symmetry φ 2 → −φ 2 , which would lead to extra stable particles and change the phenomenology of the model. In the limit of small v 3 , the only gauge-scalar mixing terms are A µ 6 ∂ µ ϕ 4 and A µ 7 ∂ µ ϕ 3 so thatφ 1 ϕ 1 andφ 2 ϕ 2 .
The mass matrix for the (pseudo)scalar fields reads: We see thatφ 3 does not mix with the other states and is a mass eigenstate. The other mass eigenstates are obtained by diagonalising the upper 3 × 3 sub-matrix. For further simplification, we will assume that the (1,2) and (2,3) entries of m 2 CP −even are much smaller than the other matrix elements, which can be achieved with sufficiently small λ H11 and λ 3 . Then ϕ 1 is approximately a mass eigenstate, which we call H (to be consistent with the notation in Ref. [25]), and m 2 H = λ 1 v 2 1 . The other two mass eigenstates are 3 with m 2 The eigenstate h 1 is identified with the 125 GeV Higgs boson and, consequently, its couplings are required to be SM-like. This translates into the requirement s θ 0.3 (see e.g. [32]).
We now turn to the vectors. In the limit v 3 v 1 , v 2 the vector sector is composed of 6 pure states which form 3 mass degenerate pairs with masses and two mixed eigenstates so that α ∈ (0 • , 60 • ). The masses are Our setup enjoys a Z 2 × Z 2 symmetry (cf. Eq. (2)). The Z 2 acts as complex conjugation, which is an outer automorphism of SU(3), while the Z 2 is a gauge transformation that acts non-trivially only on the upper entry of the SU(3) triplets. They are inherent in the Yang-Mills system and remain unbroken by interactions with matter in our minimal setting.
This discrete symmetry is in fact part of a global U (1) × Z 2 preserved by the vacuum. The global U(1) is a subgroup of the SU(3) hidden gauge symmetry and acts on the gauge fields as A µ → U A µ U † . This corresponds to (A 1,4 , A 2,5 ) → (cos βA 1,4 −sin βA 2,5 , sin βA 1,4 +cos βA 2,5 ) and leaves A 3,6,7,8 invariant. The scalar sector Eq. (5) possesses an independent global U(1) symmetry φ 1,2 → e iγ φ 1,2 . Since U acts effectively as an overall phase transformation on the scalar fields of the form Eq. (6), the vacuum preserves a combination of U(1) and U(1) . This symmetry ensures, for instance, that m A 1 = m A 2 and m A 4 = m A 5 (see [25]). Although the unbroken symmetry is U (1) × Z 2 , for our purposes it suffices to consider its subgroup Z 2 × Z 2 . The corresponding charges are given in Table 1. The lightest vector state is always A 3 . It is however not necessarily stable since |D µ φ 2 | 2 generates the coupling where allowing for the decay Hereφ 3 is produced off-shell and leads to the SM final states since the coupling ofφ 3 to h 1 , h 2 is nonzero for v 3 = 0. The masses of A 3 and χ are related by The decay . If one requires χ to be part of DM, relatively smallg necessitates therefore very small In summary, our SU (3) hidden sector adds to the particle content the following states: 8 massive vector bosons, three scalars h 2 , H,φ 3 and one pseudo-scalar χ. Given the charges under the Z 2 × Z 2 symmetry and the mass relations of Eqs. (15,18,22), different states can contribute to Dark Matter. The options are summarized in Table 2. In all cases, DM consists of 3 states. Since A 1 , A 2 are degenerate in mass, one may introduce a formal analog of the W ± bosons via the linear combinations A 1 ± iA 2 even though A 1 , A 2 have different parities. We find that such a redefinition facilitates numerical computations, in particular, what concerns the software Micromegas. This allows for the treatment of A 1 , A 2 as an effectively single (complex) DM component. A similar redefinition can be applied to another mass degenerate pair A 4 , A 5 .  (22)).
The four possible cases can be understood as follows.
is stable, because these are the lightest states with a given non-trivial Z 2 × Z 2 charge. The only possible decay would be of the type A 1 → A 2 A 3 which is kinematically forbidden. The other stable state of the hidden sector is either χ or A 3 , depending on the value of λ 4 − λ 5 . The purely vectorial DM case was studied in [25]. In this work, we will instead focus on case I, with mixed scalar-vector DM.

Multicomponent Dark Matter Phenomenology
Many of the important features of our model can be obtained by taking the limit v 1 v 2 . This reduces the number of states relevant to DM phenomenology to the DM candidates A 1,2 and χ, two mediators h 1 and h 2 , and the state A 3 whose mass is between that of A 1,2 and χ. We discuss this limit in the next subsection. Afterwards, we also consider the case v 1 v 2 where all hidden states play a role and highlight the differences between these two limits.
For v 1 v 2 , the mass scales of A 1,2 ,A 3 on one hand and A 4−7 , A 8 on the other hand are split, with the latter being higher by a factor of order v 1 /v 2 . The same happens in the scalar sector where the states H,φ 3 are parametrically heavier than h 1,2 . On the other hand, since we are interested in a relatively light χ, we take a small enough value of λ 4 − λ 5 to keep its mass below that of A 3 (cf. Eq. (22)). In practice, v 1 /v 2 3 − 5 is sufficiently large to neglect the heavier states, while we take v 1 /v 2 = 10 in our numerical studies. For , the relevant for our purposes Lagrangian is given by where we neglect the h 4 -type couplings which do not contribute significantly to the DM relic density computations. Here, the DM Lagrangian, containing the mass terms and the h 1 , h 2 interaction terms, is The couplings of h 1 and h 2 to SM matter are given by The remaining term L h-h-h represents the trilinear couplings among h 1 and h 2 , where Note that the quartic couplings λ H , λ 2 , λ H22 do not explicitly appear in the above interaction terms since they are fixed in terms of v, m h 1 , m h 2 , sin θ,g and m A : The couplings in Eq. (23) therefore are a function of the 5 new physics parameters m χ , m A , m h 2 ,g, sin θ. The hidden sector gauge couplingg acts as an overall normalization parameter for the DM interactions.
In what follows, we analyze how the DM relic density is generated as well as the constraints and prospects for Direct DM Detection.

Relic density
In conventional WIMP scenarios, the DM relic density is inversely proportional to the thermally averaged DM annihilation cross-section into SM fermions. In the case of multicomponent DM, the situation is more involved since there are additional important processes such as conversion of one DM component into another. This complicates the analysis and we therefore solve the system of Boltzmann equations numerically (cf. [17,28,33] where Y i = n i /s with n i being the corresponding number density and s being the entropy, x = m A /T and where H is the Hubble rate. The resulting evolution of the yields Y i for two benchmark parameter choices is shown in Fig. 4. In the right panel, the relic density of the two components evolves similarly to that of conventional WIMPs, i.e. it tracks the equilibrium distribution at Early times until decoupling. In the left panel, we see some modifications to this behaviour. In particular, the pseudoscalar DM components annihilates very efficiently through the h 1 resonance which depletes its energy density, while at late times the χ fraction of the DM number density increases due to the conversion process of Fig. 3. In this case, the heavier DM component gives the dominant contribution to the DM density. Our numerical analysis (see below for more details) shows that the contribution of the processes AA → A 3 h 1,2 and AA 3 → Ah 1,2 is negligible over most of the parameter space. For a qualitative discussion of our numerical results, one may thus approximate the total DM relic density by the sum of the following two contributions [37], where T f,χ , T f,A are the freeze-out temperatures of the two DM components, T 0 is the present time temperature and g eff,A,χ are the effective degrees of freedom in the Early Universe. In the second line of Eq. (32), we have used the velocity expansion σv a+2b/x (using σv a+bv 2 /3 and v 2 = 6/x, cf. e.g. [38]) and x f,i = m i /T f,i . 6 In this work we only consider the case m χ < m A as required by our model. Indeed, A 3 is always lighter than A 1,2 with our SU(3) breaking mechanism and χ must be lighter than A 3 to be stable. Hence, we include the conversion process AA → χχ, but not the reverse (at least at late times).
The relevant annihilation cross-sections are s-wave dominated, i.e. the coefficients a χ,A are not suppressed. At leading order in velocity expansion, they read • Pseudoscalar component: • Vector component: Figure 5: The ratio f A = Ω A /Ω tot in the (m χ , m A )-plane, for sin θ = 0.1, m h 2 = 500 GeV andg = 0.2 (left) respectivelyg = 1 (right). The blue, light blue, light red and red regions correspond to f A < 0.1, 0.1 < f A < 0.5, 0.5 < f A < 0.9 and f A > 0.9, respectively. In the black regions, the observed total DM relic density is correctly reproduced at the 3 σ level.
The "dark" annihilation process AA → χχ can be the most efficient A-annihilation channel since it is not suppressed by sin 2 2θ, which is subject to rather tight experimental constraints [32]. This is because the process involves only the dark sector states. As a result, the annihilation cross-section of the vector DM component is often enhanced compared to that of the scalar component.
In Fig. 5, we show the contribution of the vector component to the total DM relic density, f A = Ω A /Ω DM,tot , in the plane (m χ , m A ) with fixedg, s θ and m h 2 . We distinguish the following three regions: f A < 0.1 (blue), 0.1 < f A < 0.5 (light blue), 0.5 < f A < 0.9 (light red) and f A > 0.9 (red). The correct DM relic density is only reproduced in the black regions, so the purpose of the plot is to help understand how the composition of DM evolves as a function of parameters.
Since the total DM relic density is given approximately by Eq. (32) and x f,A ≈ x f,χ , √ḡ eff,A ≈ √ḡ eff,χ , f A mostly depends on the ratio of the pair annihilation cross-sections of the two DM components: An obvious feature of Fig. 5, which follows immediately from the above equation, is that the A µ DM component dominates when χ annihilates resonantly, and vice versa. For the regions away from the resonances, a closer inspection of σv χ σv A is required. Let us identify qualitative features of σv χ σv A . For the mass range shown in the plot, σv χ is dominated by σv χχ→bb for m W > m χ and by σv χχ→W W for m W < m χ . 7 σv A has contributions from annihilation to both dark and visible sector final states. Which one dominates depends mostly on the ratio tan θ/g. For instance, one has The ratio σv AA→bb σv AA→χχ has an additional m 2 b /m 2 χ suppression factor.
• From Eq. (36) we see that in the right plot, whereg sin θ, the dark annihilation AA → χχ dominates in most mass regions. An exception is the region where A µ is not much heavier than χ so that the dark annihilation is phase-space suppressed.
-For m χ > m W , the ratio σv χ σv A becomes For most parameter ranges of interest, the A µ annihilation cross section is much larger than that for χ. As a result, f A < 0.1. This does not however apply to the region m 2 A m 2 h 2 /4 (upper part of the plot) in which case the factor (2m A /m h 2 ) 4 can compensate the small ratio (tan θ/g) 2 .
-For m χ < m W , the χ annihilation cross section is suppressed by the b-quark mass. The resulting σv χ σv A and f A are small unless χ annihilates resonantly.
• In the left plot, sin θ andg have similar sizes so that the visible and dark A µ annihilation channels play in general comparable roles. Two features are clearly visible: resonant A µ annihilation or b-quark mass suppression of σv χ for m W > m χ lead to small f A . We find that the correct relic DM density typically requires sizableg and s θ . If g is too small, the χ DM component is overproduced. As seen in Fig. 5, atg = 0.2 one must resort to resonant χ annihilation to keep its density under control. The DM composition is very sensitive to the exact χ-mass in this case. With a larger gauge coupling,g = 1, the correct relic density is achieved in substantial regions of parameter space. We find numerically that both DM components can be as heavy as a few hundred GeV, whilegs θ 0.01 is required to keep the χ-annihilation efficient. While qualitative features of the plot can be understood semi-analytically, we have performed our numerical analysis using the software Micromegas [39] which is well suited for 2 component DM.
In Fig. 6, we show the contours of correct DM relic density in the (m A ,g)-plane (left panels) and (m χ ,g)-plane (right panels). The color coding along the contours indicates the value of f A : the red (blue) end of the spectrum refers to vector (pseudoscalar) dominance.
Many features of the plots can be understood qualitatively. In the upper left panel, the dark annihilation process AA → χχ is important, yet the resulting χ states annihilate very efficiently through the h 1 resonance into the SM fields. As a result, DM is mostly vector (apart from the small region m A m h 2 /2). The necessaryg at sin θ = 0.1 is smaller than that in [25] due to the availability of the dark annihilation channel, albeit it remains in the same ballpark of O(10 −1 ). In the right upper panel, the χ mass moves a bit further from the center of the h 1 resonance, which changes the DM composition and requires somewhat larger gauge couplings. Nevertheless, the resonance is still efficient and allows one to obtain the correct relic density with a relatively smallg.
In the lower panels, the relic density band has the resonant structure similar to that of [25]. The narrow h 1 resonance at m χ m h 1 /2 is followed by a much broader 8 resonance around m h 2 /2. The kinks in the band represent new annihilation channels becoming kinematically available, e.g. χχ → h 1 h 1 . In most regions away from the tip of the resonance, DM is predominantly pseudoscalar.
Besides the prospects for direct detection, which will be discussed in the following subsection, Fig. 6 displays the limits from perturbativity of the quartic (λ i < 4π) and gauge (g 2 i < 4π) couplings as well as those from the invisible decay of the 125 GeV Higgs boson. While the former has almost no impact on the region with the correct DM relic density, the latter excludes light values of m χ below approximately 50 GeV.

Direct detection
In this subsection, we discuss the limits from the LUX experiment [40], as well as prospects for direct detection in XENON1T [31]. The interactions of the pseudoscalar and vector DM components with nuclei are vastly different, thus it is convenient to discuss them separately.
• Scattering of the χ component: The spin-independent (SI) χ-nucleon scattering cross-section vanishes at treelevel in the limit of low momentum transfer: Thus, the pseudoscalar DM component appears to hide from detection albeit in a different manner compared to the known mechanisms which rely on the pseudo-scalar/axial-vector mediators [41,42]. The reason for it is a cancellation between the t-channel h 1 and h 2 contributions which follows from the coupling as well as the h 1 , h 2 couplings to SM matter.
Since this is an important feature of this model, let us discuss the origin of this 'blind' spot in the χ-N scattering in more detail. To this end, let us consider now the χ-nucleon interaction the interaction basis, i.e. before diagonalising the scalar mass matrix. The pseudoscalar χ interacts with the scalars ϕ 1 and ϕ 2 of the dark sector and h of the Standard Model, while only the latter couples to quarks. In the interaction basis, the effective χχN N coupling is with Here κ χ represents the χ couplings to h,ϕ 1 and ϕ 2 ; κ f gives the fermion couplings of h,ϕ 1 and ϕ 2 ; m 2 is the upper left 3×3 block of the CP-even state mass matrix given in Eq. (12). One now easily finds that We see that the reason for σ χN 0 is that we have taken λ H11 , λ 3 to be negligible. In other words, we have assumed that only one scalar mixing is significant, that is, between h and ϕ 2 , while the ϕ 1 − ϕ 2 and h − ϕ 1 ones are very small. Although this is just a simplifying assumption, it is meaningful as one does not expect all the couplings to be equally significant. The corresponding region of parameter space represents an "alignment limit" where the 3 × 3 mass matrix turns effectively into a 2 × 2 one. This yields a simple calculable model, which could perhaps be justified in a framework of a more sophisticated UV completion. Were we to relax our assumption, we would get contributions which are suppressed by the ϕ 1 − ϕ 2 and h − ϕ 1 mixing angles.
• Scattering of the A µ component: The t-channel exchange of h 1 , h 2 leads to the following SI scattering cross-section on nuclei: where µ AN = m A m N /(m A + m N ) and parametrizes the Higgs-nucleon coupling. (For an up-to-date determination of f N see e.g. [43].) In the above expression, y q are the SM Yukawa couplings, f N Tq denotes the contribution of quark q to the mass of the nucleon N and f N T G = 1 − q=u,d,s f N Tq .
Since the SI scattering of χ on nuclei is suppressed, the Direct Detection limits are obtained by comparing the experimental limits to the rescaled cross-section f A σ N A . The current limits from the LUX experiment and the projected sensitivity of XENON1T are shown in Fig. 6 by green and orange contours, respectively, assuming the exposure time considered in [40]. As discussed in the previous subsection, the χ component typically dominates the relic DM density which renders the current and future Direct Detection constraints weak to irrelevant.
Furthermore, even the regions dominated by the vector DM component are hard to probe unless one employs 1 ton detectors and a few years of exposure. This is in contrast to "typical" Higgs portal DM models (see e.g. [25]). One reason for the difference is that our setup allows for dark annihilation AA → χχ which can be dominant. The presence of this additional channel lowers the gauge couplingg required by the correct relic abundance thereby diminishing the relevance of Direct DM Detection. In addition, a low value of sin θ provides another suppression factor compared to the analysis of [25].
Finally, let us note that the unusual shape of the LUX/XENON constraints in Fig. 6 is due to the non-trivial composition of dark matter. For instance, keeping m A and m χ fixed while increasingg changes the DM composition factor f A . At large enoughg, the dark annihilation channel typically dominates which makes DM mostly pseudoscalar and thus not prone to Direct Detection. This feature is clearly visible in the plots.

Case v 1 v 2
In this subsection, we repeat our analysis for v 1 v 2 . More specifically, we take v 1 = 1.2 v 2 in our numerical studies. The main difference from the previous case is that all hidden gauge bosons have comparable masses now, cf. Eq. (15). Also the scalars H,φ 3 are expected to be as heavy as h 1 , h 2 . However, we will focus on the parameter region where H,φ 3 are heavier than the other scalars and their effect can be neglected for our purposes. This is a simplifying assumption which makes our numerical analysis tractable.

Relic density
Although the general structure of the Boltzmann equations (30) is not altered, the larger number of processes makes a semi-analytic treatment very complicated in the case v 1 ∼ v 2 . Therefore we only perform the full numerical analysis with Micromegas. Compared to the v 1 v 2 case, the following additional processes occur: • The gauge bosons A 4−7 , A 8 , which are not decoupled now, act as additional mediators of annihilation processes and therefore can enhance the annihilation rates of the vector DM component. 9 • Kinetic mixing terms give rise to additional interactions which scale approximately as m 2 χ /m 2 A . Their impact is thus limited unless the two DM components have similar masses.  Fig. 6 for v 1 = 1.2 v 2 .
• Self-interaction of the A i states could a priori lead to a sizeable effect. Our numerical analysis shows, however, that this is not the case.
In Fig. 7, we show the regions of correct DM relic density, for the same sets of parameters as in Fig. 6 (apart from v 1 = 1.2 v 2 ). We see that the isocontours of correct relic abundance do not differ substantially from those for the case v 1 v 2 .

Direct detection
As seen from Fig. 7, the Direct Detection limits change substantially. Even though there is a cancellation in σ χN as described before, for v 1 ∼ v 2 it is incomplete. The mixing term A µ 6 ∂ µ χ is now important since A µ 6 does not decouple. Eliminating this term by field redefinition leads to an additional coupling that scales as m −2 A 6 . The resulting χ-N scattering cross section is then where σ AN is given by Eq. (43). Unlike in the case v 1 v 2 (i.e. m A 6 m A ), this cross section is significant. Note that σ χN is suppressed by the factor m 2 χ /m 2 A with respect to σ AN .
In the presence of non-negligible scattering cross-sections for both DM components, the analysis of the Direct Detection limits is not straightforward. Such limits are normally given in terms of the DM-nucleon scattering cross-section as a function of the DM mass. In our case, one should compare directly the experimental outcome, i.e. the distribution of events with respect to the recoil energy, with the theoretical prediction Here with ρ 0 being the experimental value of the local DM density; m R is the reduced mass of the DM-nucleus system, F i (E R ) 2 is the form-factor due to the finite size of the nucleus (normalized to F i (0) 2 = 1), and f (v i ) is the DM velocity distribution in the detector frame. A detailed discussion of the Direct Detection limit interpretation for multicomponent DM is given in [19,20]. Here we have adopted a simple approximate procedure. We have computed the total number of recoil events, obtained by integrating the distribution of Eq. (46) 10 over a suitable range of recoil energies and multiplied the result with the number of nuclei and the exposure time in a given experiment. Given the design similarity between the LUX and XENON1T experiments, we have assumed the upper limit of 3 events for both (with two years of exposure time) [44]. This number takes into account the detector efficiency which is set to 1 in Micromegas. It is seen from Fig. 7 that the contribution from the pseudoscalar component tightens the limits from DD. Yet, the thermal DM relic density band is still out of reach of LUX. The relevant Direct Detection suppression factors include a low value of sin θ, m 2 χ /m 2 A for a light χ component as well as a relatively smallg in the domain of the broad resonance m χ ∼ O(m h 2 /2). We find that these factors are efficient enough to make the detection of a light χ beyond the reach of XENON1T, while some regions with a heavier χ can be probed. This differs from the pure vector DM case considered in [25].

Complementarity of direct and indirect detection
In this subsection, we briefly explore the possibility of observing one DM component in Indirect Detection (ID) experiments and the other one through Direct Detection.
As seen in Eqs. (33,34), the pair annihilation cross-sections are s-wave dominated and suffer no velocity suppression. Therefore, both the pseudoscalar and the vector DM components can potentially generate an ID (photon) signal from thē bb,tt, W + W − , ZZ, hh final states. However, as explained in the previous subsections, the vector component often annihilates into χχ most efficiently. As a result, the ID signal would be suppressed and thus only the pseudoscalar component could potentially be detected. This situation reverses in Direct Detection since the σ χN cross section is too small.
While a detailed analysis is beyond the scope of this work, we illustrate this point with the following example (Fig. 8). We take m h 2 = 850 GeV, m A = 450 GeV, s θ = 0.1 and focus on the range 50 − 300 GeV for the χ mass. These parameters are chosen in order to have a large pseudoscalar component since the ID rate scales with the square of the DM density.
We see that there are two small regions in Fig. 8 where both ID and DD signals could be detected. The first one is close to the h 1 resonance, i.e. for m χ 60 GeV. This region may in fact be compatible with the Galactic Center gamma-ray excess [48] although reproducing it in vicinity of the s-channel resonances is in general rather contrived [49]. The second region corresponds to masses m χ ∼ 170 − 230 GeV. In this region, the density fraction of the vectorial DM component is very low, for example, f A ≈ 0.02 at m χ = 170 GeV. This is compensated by the high DD cross section since the correct relic density requiresg 1. Specifically, for m A = 450 GeV and m χ = 170 GeV, one has f A σ N A ≈ 7.5 × 10 −47 cm 2 .
A complication here is that it is very difficult to prove that DD and ID signals come from particles with different masses. One obstacle is the large uncertainty (hundreds of GeV) in the DM mass determination through Direct Detection (cf. eg. [50]). This stems from the very weak dependence of the spectrum of recoil events on the DM mass (for heavy DM). Thus, in practice it would be challenging to prove that Dark Matter is indeed multicomponent.
Further information which can help deciphering the DM composition would be provided by collider experiments. In particular, in certain kinematic regimes, e.g. m h 2 > 2m A and m h 2 < 300 GeV, the LHC monojet events with missing energy will be able to probe the hidden sector gauge coupling in the range O(10 −1 ) − O(1) [51]. Similar constraints are obtained in Vector Boson Fusion [52] (see also [53]). Other channels can provide further probes, which will be studied elsewhere.

Conclusions
We have studied a simple UV complete set-up which entails naturally multicomponent Dark Matter with spin-1 and spin-0 constituents. The symmetry that stabilizes DM is not put in by hand, but is instead inherent in the Yang-Mills system. The model belongs to the Higgs portal category with the hidden sector consisting of SU(3) Yang-Mills fields as well as the minimal Higgs content to break this symmetry completely. Upon spontaneous gauge symmetry breaking, the system retains a global U (1) × Z 2 symmetry (assuming unbroken CP in the hidden sector). We focus on its discrete subgroup Z 2 × Z 2 which can be regarded as a DM stabilizer making the lightest vector fields and a pseudoscalar stable. These play the role of multicomponent Dark Matter.
Even though the theory is rather simple in the UV, the DM phenomenology is very rich offering a number of qualitatively different parametric regimes. For instance, the "dark annihilation" channel, where the heavier DM component pair-annihilates into the lighter component, can play an important role. Dark Matter can be mostly spin-1, mostly spin-0 or mixed. An attractive feature of the model is that the Direct DM Detection rate is suppressed as long as the SM Higgs mixes predominantly with a single scalar of the hidden sector. This phenomenon is qualitatively different from the known DD suppression mechanisms. We find that in many regions of parameter space, the Direct Detection rate is well below the LUX2016 (and sometimes XENON1T) constraint while still consistent with the thermal WIMP paradigm.
This shows, in particular, that the Higgs portal Dark Matter framework offers a number of viable options and the WIMP paradigm is not necessarily in crisis.