Singlet Scalars as Dark Matter and the Muon g-2 Anomaly

We explore a simple model containing two singlet scalars as dark matter supplemented with a vector-like lepton. Evading collider and dark matter constraints, the model is able to accommodate the correct dark matter relic abundance and the muon g-2 anomaly when the masses of the three new fields are around 100 GeV and the values of the minimum amount of required couplings are order one.

This hint of new physics beyond the Standard Model (BSM) has motivated different extensions to the SM to account for it. In particular, new interactions involving the leptonic sector along with dark matter (DM) candidates seem an appealing and exiting avenue [25][26][27][28][29][30][31][32][33][34][35][36][37][38][39][40][41]. One of the simplest extensions to the SM to account for DM and the (g − 2) µ anomaly involves the introduction of two new fields, i.e. one real scalar plus one vector-like lepton (VLL), both odd under the same discrete symmetry Z 2 (S 1 → −S 1 and ψ → −ψ). However, this set-up has been shown to * bastian.diaz@tum.de † karim1.ghorbani@gmail.com give an under-abundant relic density [33]. Extensions with three fields have also been explored [28,33], and in particular, it has been shown that extending the minimal scenario (singlet scalar plus a VLL) by including a SU (2) L scalar doublet may accommodate (g − 2) µ and DM constraints, provided that the scalar DM masses be above 500-600 GeV, and in the TeV scale when the mixed scalars are degenerated [33]. In this work we explore a simpler scenario in which we extend the minimal model (i.e. one singlet scalar S 1 and one VLL) by adding another singlet scalar S 2 , which along with S 1 both play the role of DM. The model is able to explain dark matter and (g − 2) µ at the EW scale. Multicomponent singlet scalar DM has been an interesting subject to explore [42][43][44], and in particular it has been shown to account for other event excess, such as the cosmic-ray positron-to-electron ratio [45].
The paper is organized in the following way. In section 2 we introduce the model and its freeze-out mechanism. In section 3 we present the most relevant constraints on the model. In section 4 we present different scenarios to account for (g − 2) µ and DM, and finally in section 5 a discussion with the conclusions.

Lagrangian
We introduce two real singlet scalars S 1 and S 2 , and a fermion ψ being singlet under SU (2) L with Y = −1. Two discrete symmetries are imposed Z 2 × Z 2 , where the new fields transform as In this way, the Lagrangian is given by All the parameters in the model are real. For the analysis of this paper the quartic couplings, λ i , play no role, so we can set λ i as arbitrary values in our computations. Additionally, none of the scalar fields get vev, and due to the two discrete symmetries, both S 1 and S 2 are stable. After electroweak symmetry breaking (EWSB), the masses of both singlet scalars at tree level are given by

Relic abundance
As the model contains two DM components, the total relic abundance is given by the sum of its partial contributions, i.e. Ωh 2 = Ω 1 h 2 + Ω 2 h 2 , whose value must be contrasted with the Planck measurement Ω c h 2 = 0.12 [46]. In order to achieve the correct relic abundance and at the same time respecting the new measurements on (g − 2) µ , in our framework it requires to have the following mass hierarchy The first inequality is necessary in such a way that DM conversions of the type S 1 S 1 → S 2 S 2 be effective in the thermal decoupling history, yielding Ω 2 > Ω 1 , and in this way Ω 2 covering the missing DM in the minimal scenario of two-fields scalar plus VLL [33]. The second inequality is necessary in order to have S 1 stable, otherwise it could decay into the pair ψ ± µ ∓ . Introducing x = µ/T , with µ = m 1 m 2 /(m 1 + m 2 ) and T the SM temperature, and Y i = n i /s, the coupled Boltzmann equations for the system of the two singlet scalars in thermal equilibrium in the early universe are given by with r ij (x) := Y ie (x)/Y je (x) and where σ iijj v (x) is the Maxwell-Boltzmann thermal average annihilation cross section times the relative velocity, s(T ) and H(T ) the entropy density and Hubble rate, respectively, and X stands for the SM particles, including the Higgs boson (for more details about the Boltzmann equations and its components see [47]). One important point here is that the Boltzmann equations in eq. (6) do not include coannihilation processes, e.g. S 1 + ψ ± → µ ± + V , with V = γ, Z, due to the fact that they become effective provided that m 2 ψ /m 2 1 1.2 [48]. We will see in the next section that this mass region is excluded by soft lepton searches, and therefore making coannihilation processes ineffective in the parameter space of our interest.
To illustrate the mechanism of DM production (also known as assisted-freeze-out mechanism [47,49]), in Fig. 1 we exemplify the thermal evolution of the yields of the two scalar singlets. We choose the three masses and κ 1 in the typical ball park to account for ∆a µ : m ψ = 150 GeV, m 1 = 130 GeV, m 2 = 125 GeV and κ 1 = 2.8. For simplicity we set λ H1 = λ H2 = 0. The yields are shown for λ 12 = 0.1 (dashed lines) and λ 12 = 1 (solid lines). As it can be seen in the plot, after S 1 (solid blue lines) decouples from the thermal bath (dotted lines), it continues annihilating into S 2 (red lines), also known as DM conversions, reducing its abundance even more. As a consequence, the scenario always predicts Ω 1 Ω 2 . As a matter of comparison, we show the thermal yield of a single DM scalar in a VLL portal (see Sec. 1) with the green dashed line, making explicit its sub-abundance [33]. One important aspect of this type of mechanism is the relevance of the mass shift ∆m ≡ m 1 − m 2 , which must not be bigger than a few GeV [47]. As it can be noted from the second line of eq. (6), the relevant term depending on the mass shift is This implies that the decreasing rate of Y 2 is not only regulated by λ 12 , but it also depends strongly on ∆m and µ. This is, if ∆m → 0, there is no exponential suppression, then enhancing the decreasing of Y 2 , and on the other hand, higher scalar masses (bigger µ) reduces Y 2 more ef- ficiently. This last aspect will be part of the analysis in Sec. 4. Finally, we solved the system of Boltzmann equations 6 with our own code in Python, and cross-checked our results using Micromegas code [50].

Constraints
On the scenarios in the present work we take into account the following constraints: • Dark Matter. First, we set that the total relic predicted by the model Ωh 2 must fulfill the Planck measurement Ω c h 2 = 0.12 [46], with a tolerance of ±10% accounting for uncertainties in our computation. Secondly, direct detection (DD) bounds for DM candidates with masses of O(100) GeV are set by Xenon1T [51]. Two types of portals may be present in the model, and each signal must be weighted by the relative abundance of the corresponding DM component.
DM-nucleon interactions via the VLL portal start at two-loops, but this contribution turns out to be well below Xenon1T bounds [52,53]. In the case of Higgs portal interactions the spin-independent DMnucleon scattering amplitude is given by [54,55] σ SI = λ 2 Hi f N 4π where µ = m n m i /(m n + m i ), f N = 0.3 and m n the nucleon mass. Finally, the leading indirect detection (ID) signals today are S 1 S 1 → µµ(+γ) or from the VLL and Higgs portal, respectively. We used the most recent Fermi-Lat bounds given in [56].
• Collider. LEP/LHC constrains the mass of the charged VLL to be above ∼100 GeV [57]. Furthermore, the Higgs decay into a pair of muon/antimuon may give rise to important contributions at one-loop level. Considering the measured cross section for the SM Higgs boson production, the upper limit at 95% C.L. is given by Br(h → µ + µ − ) < 4.7 · 10 −4 [58,59]. Taking into account only λ H1 coupling (equivalently for λ H2 ), the one-loop decay width is given by with F (x, y) = y x(y−1) log y + y−x x(x+y+1) log(y − x), in which x = m 2 h /m 2 ψ and y ≡ m 2 1 /m 2 ψ . Finally, searches for SUSY long-lived particles and also the CMS searches for the compressed spectra of soft leptons, set limits on the combination (m ψ , m 1 ). Here we take the results of [33], where the bounds found there apply directly to our scenario.
• Flavor Physics. We take the result in 1 as the reference deviation in the muon anomalous magnetic moment. The leading contribution to (g −2) µ in this model comes from the one-loop diagram depicted in Fig. 2, and it is given by [60] ∆a where f (t) = (2t 3 + 3t 2 − 6t 2 log t − 6t + 1)/(t − 1) 4 , which is positive for any t, and f (0) = 1.

Analysis
In the following we show how the model may account for the correct relic abundance and (g − 2) µ , evading the constraints previously described. Although the parameter space of the model is we start the analysis neglecting the Higgs portal, and afterwards we discuss the effects of it when the portal is not neglected.

Higgs portals closed
We start assuming that the Higgs portals are sub-dominant, i.e. (λ H1 , λ H2 ) 1. In this case, the model becomes very simple since only the first four parameters of 12 are relevant, and the constraints presented in the last section reduce to a few ones: relic abundance, (g − 2) µ and some collider constraints. Then, for simplicity we assume that λ H1 = λ H2 = 0.
In Fig. 3 we show some predictions of the model projected on the plane (m ψ , m 1 ) in which all the constraints are fulfilled and (g −2) µ is accounted, for κ 1 = 2.5 (plot on the top) and κ 1 = 3.5 (plot on the bottom). The white region is still allowed by the existing constraints, whereas the yellow and grey regions are already ruled out by the collider bounds. The blue region is not allowed in this model because it would make S 1 unstable (S 1 → ψ ± + µ ∓ ), and the purple one is where coannihilations become effective. The green band accounts for (g − 2) µ at 1σ, with the black line near the middle corresponds to the central value.
As it can be noted from the plot in the upper panel of  Figure 3: Results in the scenario with λ H1 = λ H2 = 0, scanning over the masses of the new fermion and the heaviest scalar S 1 , for κ 1 = 2.5 (plot on the top) and κ 1 = 3.5 (plot on the bottom). The region that can explain ∆a µ to within 1σ is colored green. The blue and red curves indicate those points that produce the correct relic abundance, with the scalar mass difference ∆m and λ 12 indicated in the legend. The blue region is not allowed by instability of m 1 , and the yellow and grey regions are LHC 13 TeV and soft leptons (s.l.) searches obtained in [33] (for details see the experimental references therein). small ∆m and λ 12 ≈ 1 enter the region allowed by the experimental constraints fulfilling (g − 2) µ . In particular, red curves with ∆m = 2 GeV enter the green zone, with tiny deviations as λ 12 changes, and as ∆m increases the model predictions deviate from the allowed region, as it can be noted for ∆m = 3 GeV (blue curves). The strong sensitivity of the model predictions to the growing of ∆m can be understood due to the fact that in order to keep the correct relic abundance at higher ∆m, it is necessary to increment the reduced mass µ to compensate the exponential suppression given in eq. (8). The sensitivity of the model predictions to λ 12 peaks at high masses where the exponential suppression in eq. (8) is less strong, therefore making the abundance Y 2 being sensitive mainly to λ 2 12 . Equivalently, the plot in the lower panel of Fig. 3 considering κ 1 = 3.5 shows that the viable mass shift picks up a little bigger value, as it is illustrated for ∆m = 5 GeV (blue curves), with an increasing range for the masses m ψ and m 1 . Note that as a big portion of the coannihilation region is excluded by LHC bounds (yellow region) and it practically does not match the (g − 2) µ region (green one) we did not show the predictions of the model in that region.
In this simple scenario, as the Higgs portal is closed, the leading contributions to ID are S 1 S 1 → µ + µ − , µ + µ − γ, γγ. However, in our scenarios S 1 is always a sub-dominant DM component, typically Ω 1 h 2 ∼ 10 −3 , then the resulting rescaled average annihilation cross section at zero velocity is suppressed by a factor of order ∼ 10 −6 . We have checked that the resulting signals are some order of magnitude below the upper bounds from Fermi-LAT.

Higgs portals opened
Now we consider the case that either one or two Higgs portals are not negligible. In this way, here we sketch what deviations occur once the Higgs portals are not negligible with respect to the previous scenario. Let us start assuming that λ H1 = λ H2 = 0. It is not possible to arbitrarily increase the values of the Higgs portal, because they not only will deviate the relic abundance found in the previous scenario, but they increase the DD signal. Without loss of generality, let us consider a point of the parameter space that fulfills (g − 2) µ and LEP/LHC constraints: (m ψ , m 1 , m 2 ) = (150, 130, 125) GeV, κ 1 = 2.8 and λ 12 = 1. In Fig. 4(upper panel) we show the deviation in the relic abundances as the Higgs portal couplings change according to λ H1 = λ H2 . For λ H2 0.01, mod- In the plot above we consider λ H1 = λ H2 , whereas in below λ H2 = 0. The blue and red curves are allowed by DD bounds, whereas the grey one is the already ruled out by Xenon1T. The green region is forbidden by Br(h → µ + µ − ) measurements, and the dotted horizontal line is the measured relic abundance Ω c h 2 .
ification on Ω 2 starts appearing with respect to the case λ H1 = λ H2 = 0 and almost simultaneously the model enters the region in conflict with the DD bounds (grey curve). Therefore, provided that the Higgs portal couplings take values lower than O(0.01), no significant deviations with respect to what we have found in Sec. 4.1 will appear in this case.
Since the scalar S 2 gives most of the DM relic abun-dance (i.e. in this scenario we always have Ω 1 Ω 2 ) and therefore most sensible to the DD bounds via the Higgs portal, another interesting option is to set λ H1 = 0 keeping λ H2 very small. In this case, higher values for λ H1 can be reached before the DD bounds start to be relevant. In particular, for the same parameter space point described in the previous paragraph, in Fig. 4(lower plot) λ H1 0.1 still fulfills the relic abundance while evading DD bounds. Therefore, allowing the Higgs portals coupling take nonvanishing values, DD bounds kill all the significant deviations to the predictions obtained in Sec. 4.1 (Higgs portal closed).
As a matter of complement, note that the scenario λ H1 = 0 and λ H2 = 0 is not unstable under radiative corrections due to the following. At one-loop, the effective S 2 S 2 h coupling in the limit in which λ H2 1 is given by [47] λ eff H2 = λ H2 + λ H1 λ 12 16π 2 .
As we have seen, the maximum value of λ H1 allowed by the DD bounds is O(0.1), then implying that at one-loop λ eff H2 ≈ 6.3 × 10 −4 , for λ 12 ∼ O(1). As we have seen in this section, such small values for λ eff H2 produce no deviation on the relic abundance and evade easily the DD bounds.

Discussion and Conclusions
In this work we have shown that a simple unexplored extension to the SM containing two real singlet scalars and one vector-like lepton not only accounts for the measured DM relic abundance evading stringent experimental constraints, but it can also fulfill the recent 4.2σ deviation in (g − 2) µ , to be tested in the forthcoming years with more data. Provided that the mass of the new three states be around the Fermi scale and the mass shift between the singlet scalars is not more than a few GeV, the (g − 2) µ anomaly is perfectly achieved, with the minimal requirement that the new couplings be of order one. We have seen that this set-up implies that the inert scalar S 2 captures most of the relic abundance, whereas the sub-leading component S 1 plays the role of enhancing the anomalous magnetic moment of the muon at one-loop. Furthermore, we have seen that the Higgs portal couplings play a secondary role in this framework, although they can still take relatively sizable values in agreement with all the experimental constraints. Finally, an additional comment about the possibility that also S 2 couples to the muon via the VLL such that both scalars participate in the (g − 2) µ (this can be done assuming ψ → −ψ under both Z 2 and Z 2 ). In that case, one can have κ 1 and κ 2 with slightly lower values than the obtained ones in this work. That scenario is dangerous because the scalar with the highest DM relic abundance (i.e. at least half of Ω c h 2 ) will present strong indirect signals (e.g. bremsstrahlung), then entering in conflict with the present ID bounds.