Minimal Two-component Scalar Doublet Dark Matter with Radiative Neutrino Mass

We propose a minimal extension of the Standard Model to accommodate two-component dark matter (DM) and light neutrino mass. The symmetry of the Standard Model is enhanced by an unbroken $\mathbb{Z}_2 \times \mathbb{Z}'_2$ such that being odd under each $\mathbb{Z}_2$, there exists one right handed neutrino and one inert scalar doublet. Therefore, each of the $\mathbb{Z}_2$ sectors contribute to ($i$) light neutrino masses radiatively similar to the scotogenic models while ($ii$) the two neutral CP even scalars present in two additional inert doublets play the role of dark matters. Focussing on the intermediate range of inert scalar doublet DM scenario: $M_W \leq M_{\rm DM} \lesssim 500 \; {\rm GeV}$, where one scalar doublet DM can not satisfy correct relic, we show that this entire range becomes allowed within this two-component scalar doublet DM, thanks to the inter-conversion between the two DM candidates in presence of neutrino Yukawa couplings with dark sector.


I. INTRODUCTION
There have been irrefutable amount of evidences in favour of the existence of nonluminous, non-baryonic form of matter in the universe, popularly known as dark matter (DM). The presence of this form of matter has also been supported by astrophysical observations like the ones related to galaxy clusters by Fritz Zwicky [1] back in 1933, observations of galaxy rotation curves in 1970's [2], the more recent observation of the bullet cluster by Chandra observatory [3] along with several galaxy survey experiments which map the distribution of such matter based on their gravitational lensing effects. There is equally robust evidence from cosmology as well, suggesting that around 26% of the present universe's energy density is in the form of dark matter. In terms of density parameter Ω DM and h = Hubble Parameter/(100 km s −1 Mpc −1 ), the present DM abundance is conventionally reported as [4]: Ω DM h 2 = 0.120 ± 0.001 at 68% CL.
In spite of these astrophysical and cosmology based evidences, there have been no detection of particle DM at any experiments. The direct detection experiments like LUX [5], PandaX-II [6,7] and Xenon1T [8,9] have continued to produce null results so far. As the Standard Model (SM) of particle physics can not accommodate such a form of matter, several beyond the Standard Model (BSM) proposals have been put forward [10] out of which the weakly interacting massive particle (WIMP) paradigm is the most widely studied one.
While such interactions are capable of explaining correct relic density of DM, the same interactions can also give rise to production of DM particles at the large hadron collider (LHC) [11]. However, nothing is found so far in these searches also, putting strict bounds on DM coupling to the SM particles, particularly quarks. Another detection prospect lies in the indirect detection frontier where searches are going on to find excess of antimatter, gamma rays or neutrinos, originating perhaps from dark matter annihilations (for stable DM) or decay (for long lived DM). While no convincing DM signal has been observed yet, there are tight constraints on DM annihilations into SM particles [12], specially the charged ones which can finally lead to excess of gamma rays for WIMP type DM.
Though the null results mentioned above have not ruled out all the parameter space for a single particle DM models yet, it may be suggestive of presence of a much richer DM sector. The idea may be natural given the fact that the visible sector is made up of several generations for single type of particles. There have been several proposals for multi-component WIMP dark matter during last few years, some of which can be found in .
Such multi-component DM scenarios, even if both the DM candidates are of the same type, can have very interesting signatures at direct as well as indirect detection experiments, as studied in [42][43][44][45][46][47][48][49][50][51][52][53][54][55][56][57][58]. Since direct and indirect detection (considering annihilations only, for stable DM) rates of DM are directly proportional to the DM density and DM density squared respectively and thereby producing tight constraints on single DM models, multicomponent DM models can remain safe from being ruled out by such direct and indirect search constraints if the relative densities of different DM components are within appropriate limits. On top of that, such multi-component DM often comes with additional features like giving rise to interesting indirect detection signatures like monochromatic X-ray or gamma ray lines, as explored in several works, see for example [17,24,29,[54][55][56] and references therein.
Multi-component DM scenarios may also be connected with other sector of particle physics. One such immediate possibility evolves through a probable connection with neutrino physics, particularly with the origin of neutrino mass and mixing. Results of several experiments in last two decades like T2K [59], Double Chooz [60], Daya Bay [61], RENO [62] and MINOS [63] have confirmed the existence of non-zero but tiny neutrino mass and large (compared to quark mixing) leptonic mixing [59][60][61][62][63][64][65][66][67][68]. Similar to the case of DM, these experimental observations provide clear indication for BSM physics as neutrino mass can not be explained within SM framework. Several BSM models attempt to explain tiny neutrino mass by incorporating additional fields. Apart from the conventional type I seesaw [69][70][71][72], there exist other variants of seesaw mechanisms also, namely, type II seesaw [73][74][75][76][77], type III seesaw [78] and so on.
It is particularly interesting to think of possible connection between the origin of neutrino mass and dark matter [79,80], and perhaps such a connection is most straightforward in scotogenic scenarios, originally proposed by Ma [81]. In scotogenic type of model, the Z 2 odd particles take part in radiative generation of light neutrino masses, while the lightest  [85][86][87][88][89][90][91][92][93][94][95][96][97][98][99][100][101][102][103]. We find that the total relic abundance of two component scalar doublet DM can be satisfied in the intermediate region while being consistent with neutrino oscillation data. The model also predicts the lightest neutrino mass to be zero. While the DM candidates satisfy the constraints from direct and indirect detection experiments, there lies the tantalising possibility to probe these scenarios at future searches in these frontiers and also in rare decay experiments like µ → eγ. This paper is organised as follows. In section II, we discuss the model and the particle spectrum. In section III we discuss the details of two component dark matter pointing out the different annihilation channels contributing to the individual and total relic abundance, constraints from direct, indirect search followed by discussions on the constraints from neu- 1 In another recent work [39], two component fermion DM was proposed as a new anomaly free gauged B − L model. 2 Fermion singlet DM in such models typically require large Yukawa couplings to satisfy correct relic and often run into the problems of vacuum stability [84]. trino oscillation data in section IV. We then briefly comment on lepton flavour violation in section V. We discuss our results in section VI and finally conclude in section VII.

II. THE MODEL
We have extended the particle content of the Standard Model by introducing two SU (2) L scalar doublets η 1 and η 2 and two right handed (RH) neutrinos, N 1,2 . Furthermore, we include additional discrete symmetries, Z 2 ×Z 2 under which all SM fields transform trivially.
The charge assignments of these additional fields under SM gauge symmetry as well as additional global discrete symmetries are indicated in Table I where α, β = e, µ, τ stand for different generations of SM leptons. Note that as each RH neutrinos are odd under two different Z 2 sectors, the corresponding RH neutrino mass matrix remains diagonal.
This can be written as follows. where and In order to keep Z 2 × Z 2 unbroken so as to guarantee the stability of DM components, the neutral components of η 1 and η 2 are chosen not to acquire any non-zero vacuum expectation value (VEV) and hence they can be identified as two inert Higgs doublets (IHD). These IHDs can be parametrised as On the other hand, the neutral component of the SM Higgs field acquires a non-zero VEV (denoted by v) which is responsible for electroweak symmetry breaking (EWSB). We parametrise the Higgs field H as The scalar potential should be bounded from below in order to make the electroweak vacuum stable. This poses some constraints on the scalar couplings of the model. In addition to this, all the relevant couplings should also maintain the perturbativity. These bounds together with the unitarity limit have been studied for the three Higgs doublets and can be found in [104], from which we obtain the limits in our case.
We are effectively left with two-component inert DM scenario where the presence of light neutrinos are also taken cared of. As specified before, we consider the CP even neutral components H 1 , H 2 as two DM candidates without any loss of generality. Analogous to the single IHD scenario, we define λ L1 ≡ λ 3 +λ 4 +λ 5 2 and λ L2 ≡λ 3 +λ 4 +λ 5

2
, which denote the individual Higgs portal couplings of two DM candidates respectively. For our analysis purpose, we first implement this model in LanHEP [105] choosing the independent parameters in the scalar sector as For simplicity, we will consider the couplings of quartic interaction between two IHDs (η 1 and η 2 ) to be the same and identify it by λ 12 , i.e.λ 3 , λ 4 , λ 5 = λ 12 . We express other couplings of the scalar sector in terms of these independent parameters as follows,

III. DARK MATTER PHENOMENOLOGY
The set-up contains two dark matter components: H 1 and H 2 . In order to find the final relic density, their annihilations and co-annihilations are to be considered. In addition, the role of neutrino Yukawa couplings are also important. Below we provide a systematic approach to calculate the relic density in our scenario. Constraints from dark matter search are also mentioned.

A. Relic Density
For a single component WIMP type DM, the DM candidate with mass m DM is part of the thermal plasma in the early universe which eventually freezes out when the rate of annihilations fall below the rate of expansion of the universe. The final abundance can then be obtained by solving the Boltzmann equation for the DM number density. In fact, for DM annihilations dominated by s-wave processes only, the relic abundance can be approximated as [106] Ω DM h 2 = 1.07 × 10 9 GeV −1 where g * and g * s are the effective relativistic degrees of freedom that contribute to the energy density and entropy density, respectively. x f is to be determined from the parameter x = m DM /T evaluated at the freeze-out temperature T f and is given by with g being the number of internal degrees of freedom of the DM. The thermally averaged annihilation cross section, given by [107] is evaluated at T f and denoted by σv f . The freeze-out temperature T f is derived from the equality condition of DM interaction rate Γ = n DM σv with the rate of expansion of the universe H(T ) In the above expression of Eq.(10), K i (x)'s are the modified Bessel functions of order i. As is well known, if the mass splitting within the DM multiplet is relatively small, there can be additional contributions from coannihilation channels [108], whose importance in IDM has already been discussed in several earlier works.
In the present model, we have two dark matter candidates H 1 and H 2 . Since both the candidates contribute to the dark matter relic density obtained from Planck [4] experiment, one must satisfy the following relation: where h denotes the reduced Hubble parameter and the relic density of the H 1 and H 2 is given by Ω 1 h 2 and Ω 2 h 2 respectively. Since there exists annihilation channels through which

FIG. 2. Coannihilation Channels
Here y i (i = H 1,2 ) is related to Y i by y i = 0.264M P l √ g * µY i and similarly for equilibrium where the equilibrium distributions are now recasted in terms of µ having the form 3 We adopt the notation from a recent article on two component DM [31] Here M Pl = 1.22 × 10 19 GeV and g * = 106.7 and X represents SM particles. One should note that the contribution to the Boltzmann equations coming from the DM-DM conversion (corresponding to Fig. 3) will depend on the mass hierarchy of DM particles. This is described by the use of Θ function in the above equations. These coupled equations can be solved numerically to find the asymptotic abundance of the DM particles, y i which can be further used to calculate the relic: where x ∞ indicates a very large value of x after decoupling.
In presence of neutrino Yukawa couplings, the Boltzmann equations get modified and are given by The total relic density of DM follows from the combined contribution of both the components and is given by As mentioned earlier, DM parameter space can be constrained significantly by the null results at different direct detection experiments such as LUX [5], PandaX-II [6,7] and Xenon1T [8,9]. There are two ways scalar DM can scatter off nuclei at tree level in our model.
One is elastic scattering mediated by SM Higgs boson while the other one is the inelastic one mediated by electroweak gauge bosons. The latter can be kinematically forbidden by considering large mass splitting between IHD components, typically larger than the average kinetic energy of DM particle. The spin independent elastic scattering cross section mediated by SM Higgs (shown in Fig. 5) is given as [88] σ where µ i,n = m n m i /(m n + m i ) is the DM-nucleon reduced mass and λ L i is the quartic coupling involved in DM-Higgs interaction. The index i stands for DM candidate in our scenario: H 1 , H 2 . A recent estimate of the Higgs-nucleon coupling f gives f = 0.32 [109] although the full range of allowed values is f = 0.26 − 0.63 [110]. In this two-component DM framework, the spin-independent cross section relevant for each of the candidate can be expressed as Latest results from Xenon-1T experiment provides a strong constraint on single component IDM as it restricts λ L coupling significantly. However due to the presence of two DM components here in our set-up, such tight constraints can be evaded by suitable adjustment of relative DM abundance. We will discuss the status of our model at direct detection frontier in subsequent sections.

C. Indirect Detection
As mentioned earlier, WIMP DM candidates have good prospects at indirect detection experiments looking for excess of gamma rays. DM particles can annihilate and produce SM particles, out of which photons (and also neutrinos), being electromagnetically neutral, have better chances of reaching the detector from source without getting deflected. Following the notations of [12], the observed differential gamma ray flux produced due to the DM annihilation can be computed as where Ω is the solid angle corresponding to the observed region of the sky, σv is the thermally averaged DM annihilation cross section, dN/dE is the average gamma ray spectrum per annihilation process and the astrophysical J factor is given by In the above expression, ρ is the DM density and LOS corresponds to line of sight. Therefore, measuring the gamma ray flux and using the standard astrophysical inputs, one can constrain the DM annihilation first into different charged final states like µ + µ − , τ + τ − , W + W − , bb which in turn produces the gamma rays. As discussed in case of direct detection, here also our set-up carries some flexibility as far as indirect detection constraints are concerned.
We incorporate the global analysis of the Fermi-LAT and MAGIC observations of dSphs [12] for this purpose. The bounds quoted in [12] consider 100% annihilation of DM into particular final states as well as assume a single DM component which fills the entire 26% of the universe. Since our construction involves deviation from these consideration, we can make the bounds weaker by playing with the branching fraction to a particular final states and simultaneously changing the relative fractional abundance. This is because the DM annihilation rates are directly proportional to number density squared of DM in local neighbourhood.
The light neutrino mass matrix m ν can be diagonalised through where m diag where c ij = cos θ ij , s ij = sin θ ij and δ is the leptonic Dirac CP phase. The diagonal phase matrix U P = diag(1, e iα , e iβ ) contains the Majorana CP phases α, β which remain undetermined at neutrino oscillation experiments. We summarise the 3σ global fit values in Table II from the recent global fit [111], which we use in our subsequent analysis. Note that the two Z 2 sectors in our model can generate at most two light neutrino masses 4 , which needs to be consistent with light neutrino data mentioned above. This leaves us with two possibilities (a) m 1 = 0, m 2 < m 3 , (b) m 1 < m 2 , m 3 = 0 corresponding to normal and inverted hierarchies respectively. Since the inputs from neutrino data are only in terms of the mass squared differences and mixing angles, it would be useful for our purpose to express the Yukawa couplings in terms of light neutrino parameters. This is possible through the Casas-Ibarra (CI) parametrisation [113] extended to radiative seesaw model [114] which allows us to write the Yukawa couplings as where Λ is the 2×2 diagonal matrix with eigenvalues defined in Eq. (22) and R is a complex orthogonal matrix [112,115,116] having the form, Using Eq.(25), the elements of 3 × 2 Yukawa matrix can be obtained with specific choices of the complex angle z. The same calculation can be repeated for inverted hierarchy as well.
In the subsequent sections, we discuss how the constraints from neutrino sector can play a non-trivial role in the dark matter parameter space from relic abundance criteria. 4 With the involvement of a third RH neutrino to start with and making it very heavy, it effectively leads to the two RH neutrino scenario we consider here. As a consequence of this limit, it is shown [112] that one of the light neutrino remains massless.

V. LEPTON FLAVOUR VIOLATION
Since the charged lepton flavour violating (LFV) decays are very much suppressed in the SM, any observation of such effects will be a clear signature of beyond the SM physics. In our model, due to the coupling of each Z 2 sector particles (N i and η i ) to the SM leptons, one may expect some contribution to such LFV effects at one loop level. The same fields that take part in the one-loop generation of light neutrino mass as shown in Fig. 6 also mediate LFV processes like µ → eγ. The neutral scalar in the internal lines of Fig. 6 will be replaced by their charged counterparts (which emit a photon) whereas the external fermion legs can be replaced by µ, e respectively, generating the one-loop contribution to µ → eγ.
As the couplings and masses involved in this process are the same as the ones that generate light neutrino masses and play a role in DM relic abundance, we can no longer choose them arbitrarily. It should be noted that each Z 2 sector contributes separately to this process and hence we have to add the respective contributions at amplitude level. Adopting the notations from [117], we can write Here α em = e 2 /4π is the electromagnetic fine structure constant, G F is the Fermi constant, MEG experiment provides the most stringent upper limit on the branching ratio: Br(µ → eγ) < 4.2 × 10 −13 [118].

VI. RESULTS AND DISCUSSION
As mentioned earlier, we write the model in LanHEP     We also check the effects of considering different mass splittings ∆M 1 = ∆M 2 . As we increase ∆M 2 to 10 GeV while keeping ∆M 1 at 1 GeV, the relic abundance of H 2 decreases, a feature observed in single component IDM also in the high mass regime. The results are shown in Fig. 9 (considering m H 2 > m H 1 ) which shows that the masses of H 2 are shifted to higher mass regions in order to satisfy the total DM relic abundance. This is expected as we know based on our knowledge of single component inert doublet DM analysis that making the mass splitting related to one IHD more, the corresponding relic density would be less. This clarifies why there is a separation between the green and orange patches (with λ 12 = 0). Similar situation prevails when nonzero λ 12 value is switched on as well, which can be seen from black (brown) and purple (red) patches. Also notice that the existence of the symmetric point similar to the case with equal ∆M i (about which the relic contour lines gain a typical shape with non-zero λ 12 ) is lost here. This is related to the fact that now two components of DM have different type of co-annihilations and hence contribute to the total relic differently. This is also reflected from a mixed up distribution pattern of correct relic satisfied points put in m H 1 , m H 2 plane in the right panel of Fig. 9. One more point to notice is that there is a shift toward the higher mass range of DM as compared to the case displayed in Fig. 7 with m H 2 > m H 1 .  LZ [121], XENONnT [122], DARWIN [123] and PandaX-30T [124].
We also check the prospects for indirect detection of DM in our model by specifically focussing on W + W − final states from DM annihilations. This is due to the chosen mass range of DM where annihilation to this final states is the most dominant one. We incorporate the relative abundance of the two DM candidates while calculating their annihilation rates.
We include the factor coming from the branching fraction to this final state from DM annihilations. Our findings are displayed in the right panel of Fig. 10 This has non-trivial connections to the light neutrino sector as the same Yukawa couplings also play a role in generating light neutrino masses. The new annihilation channels that will come into play now are the ones shown in Fig. 4. We now have an additional free parameter which is important for dark matter phenomenology defined as ∆M N iHi = M i − m H i . For analysis purpose, we fix it to ∆M N iHi = 10 GeV value. Such small mass splitting between neutral singlet fermion and DM candidates enhances the co annihilation cross section effectively having significant effect on the relic abundance as we will see soon.
In obtaining the relic density and direct/indirect detection cross sections in our model,   However while making both of these splittings smaller, it will increase branching ratio related to the prediction of LFV decays as seen from the Fig. 13. This is particularly due to the increase in Yukawa couplings for smaller values of mass splitting, in accordance with the CI parametrisation. Hence with not-so-small values of the mass splittings, we can conclude that the range of the parameter space allowed by the DM relic density, direct and indirect search results and neutrino data is mostly unaffected by the LFV constraints. The two dark matter candidates are chosen to be the neutral components of the two scalar doublets which are thermally produced in the early universe by virtue of their electroweak gauge interactions. We particularly focus on the mass range in between W boson mass and approximately 500 GeV where single component scalar doublet dark matter can not satisfy correct relic abundance. While scalar doublet relic does not depend crucially upon the Yukawa couplings to leptons, turning the Yukawa interactions on has interesting consequences for relic due to additional coannihilation channels. The inter-conversion between two dark matter candidates also play instrumental role in deciding the total relic abundance. We show that the model can satisfy correct total relic abundance criteria in the intermediate mass regime in agreement with the latest data from Planck mission. In addition to this, the model remains very much predictive at ongoing direct, indirect detection experiments as well as rare decay experiments looking for charged lepton flavour violation.
Unlike typical multi-component dark matter models, our model gets more restrictive due to electroweak gauge interactions as well as the non-trivial roles both the dark matter particles play in generating light neutrino masses, opening up detection prospects through lepton portals like charged lepton flavour violation. While the origin of the discrete symmetries Z 2 × Z 2 remains unexplained in the current work, we leave a more detailed study looking for UV completion of such models to an upcoming work.