Limiting the Heavy-quark and Gluon-philic Real Dark Matter

We investigate the phenomenological viability of real spin half, zero and one dark matter candidates, which interact predominantly with third generation heavy quarks and gluons via the twenty-eight gauge invariant higher dimensional effective operators. The corresponding Wilson coefficients are constrained one at a time from the relic density $\Omega^{\rm DM} h^2$ $\approx$ 0.1198. Their contributions to the thermal averaged annihilation cross-sections are shown to be consistent with the FermiLAT and H.E.S.S. experiments' projected upper bound on the annihilation cross-section in the $b\,\bar b$ mode. The tree-level gluon-philic and one-loop induced heavy-quark-philic DM-nucleon direct-detection cross-sections are analysed. The non-observation of any excess over expected background in the case of recoiled Xe-nucleus events for spin-independent DM-nucleus scattering in XENON-1T sets the upper limits on the eighteen Wilson coefficients. Our analysis validates the real DM candidates for the large range of accessible mass spectrum below 2 TeV for all but one interaction induced by the said operators.


I. INTRODUCTION
Regardless of the unequivocal astrophysical evidences from rotation velocity curves [1] via mass-to-luminosity ratio, Bullet Cluster [2], and precision measurements of the Cosmic Microwave Background (CMB) from WMAP [3] and PLANCK [4] satellites etc., we have no clue about the fundamental nature of dark-matter (DM), which forms the dominant matter component of the Universe. Dedicated experiments for the direct-detection of DM such as LUX [5], XENON-1T [6], DarkSide50 [7], PandaX-4T [8] and CRESST-III [9], PICO-60 [10], PICASSO [11] are designed to measure the momentum of the recoiled atom and/ or nucleus due to the scattering of DM particles off the sub-atomic constituents of the detector material [12][13][14]. We haven't seen any significant signal excess over the expected background yet.
Several experiments, past and ongoing, have constrained a plethora of viable UV complete dark matter models formulated by writing renormalizable Lagrangians with heavy non-SM mediators (spin 0 ± , 1/2, 1, 2) facilitating interactions between DM and SM particles [21][22][23][24][25][26][27][28][29][30][31]. The models with Higgs bosons as mediators have been excluded by the recent collider experiments [32]. Analogous analysis has also been performed in the domain of the effective field theory (EFT) for the EW-Boson-philic DM operators [33][34][35]. Recently, the GAMBIT collaboration performed a global analysis for signatures of SM gauge singlet Dirac fermion DM at LHC in the EFT set-up with simultaneous activation of fourteen effective operators constructed in association with light-quarks, gluons and photons (up-to mass dimension seven) where the authors have disfavored the exclusive contribution of DM ≤ 100 GeV [36].
The interactions of heavy-quark-philic DM are also constrained by the ongoing collider experiments. The CMS [59,60] and ATLAS [61] collaborations investigated the dominant scalar-mediated production of t-quark-philic DM particles in association with a single top quark or a pair of tt. CMS collaboration recently conducted an exclusive search analysis for the b-philic DM in reference [62]. Many authors in the literature have explored mono-jet + / E T mode at the LHC for the viable signatures of the scalar current induced heavyquark-philic DM interactions [51,[63][64][65]. A comprehensive search analysis for heavy scalar mediated top-philic DM models for mono-jets+ / E T , mono-Z, mono-h and tt pair productions at LHC can be found in [66,67].
The phenomenology of spin 1/2, 0 and 1 real DM deserves a special mention in the context of their contributions to the DM-nucleon scattering events in the direct-detection experiments. In the absence of the contributions from vanishing vector operators for Majorana and real vector DM, the spin-independent DM-nucleon scattering cross-sections are found to be dominated by the respective scalar and dimension-8 twist operators [68][69][70][71][72][73][74][75]. The pseudo-scalar current contribution to the DM -nucleon scattering cross-section is found to be spin-dependent and velocity suppressed. Similar studies have also been undertaken for lepto-philic real DM candidates [39,40].
In this context, it is worthwhile to investigate the WIMP DM phenomenology induced by heavy-quark-philic and gluon-philic scalar, pseudo-scalar, axial-vector, and twist-2 real-DM operators and explore whether they satisfy the relic density criteria and other experimental constraints for a DM mass between 10 GeV and 2 TeV. We introduce the effective Lagrangian for real spin 1/2, 0 and 1 DM particles in section II. In section III, we discuss the phenomenology of the dark-matter. We investigate the cosmological constraints on the DM in sub-section III A. In sub-section III B, we predict and analyse the thermally averaged DM pair annihilation cross-section. Sub-section III C details the computation of DM-nucleon scattering cross-sections and expected recoil nucleus event(s) for XENON-1T [6] set up.
Section IV summarises our study and observations.

II. HEAVY-QUARK AND GLUON-PHILIC DM EFFECTIVE OPERATORS
It is a fact that in a beyond-standard model (BSM) renormalizable gauge theory, as long as the energy range of DM-SM interaction is much below the mediator mass, the study of DM interactions can be restricted by the DM and SM degrees of freedom and their symmetries.
The higher dimensional effective operators are obtained from the BSM Lagrangian by writing the operator product expansion (OPE) of currents in the limit p 2 / m 2 Med. << 1, where p µ represents the four momentum of the virtual mediator of mass m Med. . This facilitates the effective contact interaction between DM and any SM third generation heavy quark/ gluon, assuming that the mediator mass scale m Med. is of the order of the effective theory cut-off (∼ Λ eff ), which is much heavier than the masses of the SM and DM fields in general.
For example, in renormalizable models, the interaction between a Majorana DM χ and a third generation quark ψ can be written as (1) where η and ζ µ are electrically charged scalar and vector fields respectively. This interaction Lagrangian facilitates the Majorana DM -quark interaction via t-channel exchange of the heavy η and/ or ζ µ . Thus, expanding the propagator in powers of p 2 /m 2 Med , yields the higher-dimensional effective four-fermion interaction for Majorana DM and heavy-quarks.
The scalar, pseudo-scalar, axial-vector, and twist-2 operators are all induced by this expansion [68]. When compared to a spin-2 graviton mediated model, the twist operator corresponds to the contribution from traceless part of the energy-momentum tensor T µν [76]. In various WIMP-inspired renormalizable electroweak models, such effective interactions between the SM and Majorana DM particles are also realised by s-channel processes, where the interactions are mediated by non-SM heavy scalar/pseudo-scalar/axial-vector or spin-2 tensor particles [77,78].
Since the DM-gluon interactions can be naturally realised via one-loop interactions of the DM particles either with SM heavy quarks or BSM non-singlet coloured spin 1/2, 0 and 1 exotics at the next order in the strong coupling constant, it becomes all the more necessary to include the study of effective operators constructed independently with a Majorana DM bilinear and a pair of gluons at the leading order of ∼ α s / π.
The phenomenological effective Lagrangian for the heavy-quark-philic and gluon-philic Majorana DM, χ, is written as: representing the third generation quark-philic scalar, pseudoscalar, axial-vector and twist-2 type-1 and type 2 operators respectively along with O g representing the gluon-philic scalar, pseudo-scalar and twist-2 type-1 and type-2 operators respectively are defined as The second rank twist tensor currents O q µν and O g µν for the heavy-quarks and gluons respectively are given as where D µ L and D µ R are the covariant derivatives for left and right handed quarks respectively in SM. The contribution from the vector operator vanishes for the real particles. We exclude the contribution of the dimension-9 twist-2 type-2 operators O g in our analysis because we are only interested in the effective operators up to mass dimension-8.
We extend the domain of our analysis to include the effective contact interactions of real scalar φ 0 and vector V 0 µ DM candidates with SM third generation heavy quarks and gluons. The effective scalar DM Lagrangian are given as where The O are the scalar and second rank twist-2 type-2 operators respectively.
There are no contributions from pseudo-scalar and axial vector currents for the scalar DM operators.
The effective vector DM Lagrangian is given as The heavy-quark- representing the third generation quark-philic scalar, pseudo-scalar, axial-vector and twist-2 type-2 operators respectively and gluon-philic operators corresponding to the scalar, pseudo-scalar and second rank twist-2 type-2 interactions respectively are defined as where µνρσ is the totally antisymmetric tensor with 0123 = +1.
In this study, each operator's phenomenology is investigated independently, assuming that the unique interaction contributes to the total relic density for DM.

III. PHENOMENOLOGICAL CONSTRAINTS ON REAL DM
A. Contribution of DM to the relic density The present day relic abundance of the DM species n DM (t) can be calculated by solving the Boltzmann equation where n eq DM is the DM number density at thermal equilibrium, H 0 is the Hubble constant, | v DM | is relative velocity of the DM pair and σ ann | v DM | is the thermal average of the annihilation cross-section [79]. It is customary to parametrise ρ DM ≡ Ω DM h 2 ρ c , where ρ c ≡ 1.05373 × 10 −5 h 2 /c 2 GeV cm −3 is the critical density of the universe and dimensionless h is the current Hubble constant in units of 100 km/s/Mpc. Solving the Boltzmann equation [80] we then get where parameter x F ≡ m DM /T F is a function of degrees of freedom of the DM g, the effective massless degrees-of-freedom g eff (∼ 106.75 and 86.75 above m t and m b mass-thresholds respectively) at freeze-out temperature T F : The predicted DM relic density Ω DM h 2 of 0.1138 ± 0.0045 and 0.1198 ± 0.0012 by WMAP [3] and Planck [4] collaborations respectively restrict the thermal averaged DM annihilation cross-section σ | v DM | ≥ 2.2 × 10 −26 cm 3 s −1 as a smaller thermally averaged cross-section would render large DM abundance which will over-close the universe. We have analytically calculated the DM pair annihilation cross-sections to a pair of third generation heavy quarks and gluons in the Appendix A.
We investigate the DM relic density contributions propelled by the thermally averaged annihilation cross-sections given in the appendix B corresponding to the scalar, pseudoscalar, axial-vector and twist-2 currents of the heavy quarks for a given Majorana/ scalar/ vector DM. The annihilation channels induced by the scalar, pseudo-scalar and twist-2 currents of gluons also contribute to the relic abundance of Majorana, scalar and vector DM.
Throughout our investigation, we assumed that the operators would be actuated one at a time. As a result, the constraints on all twenty-eight Wilson coefficients C q, g DM i O j / Λ n in TeV −n (DM i ≡ χ, φ 0 , V 0 and O j ≡ S/ PS/ AV/ T 1 / T 2 are derived by triggering the lone contribution from the associated operator to fulfill the relic density Ω DM h 2 ≈ 0.1198±0.0012 [4]. Due to its only action, these constraints can be regarded as a cosmologically acceptable lower-limits on the corresponding Wilson Coefficients, which translate to upper bounds on the respective cutoffs C q, g DM i O j −1/n Λ in units of TeV for a given DM mass. When the Wilson coefficient is greater than its lower-limit for a given DM mass, it partially meets the relic density.
It is important to note that because we can raise the annihilation cross-section by turning on more Wilson coefficients, we can easily imagine a scenario in which they all reside below their lower-limit without over-closing the universe with DM. To put it another way, if we set two distinct Wilson coefficients to their lower limits for a given DM mass, we plainly underproduce DM because the overall annihilation cross-section is larger. Therefore, in order to be consistent with a multiple species of DM model satisfying the relic density constraint,

FIG. 2:
We depict relic density contours satisfying Ω φ 0 h 2 = 0.1198 [4] and shaded cosmologically allowed regions in the plane defined by Majorana DM mass m φ 0 and |C The b-quark-philic, t-quark-philic and gluon-philic scalar DM contours in figures 2a, 2b and 2c are drawn corresponding to scalar and twist-2 type-2 operators respectively. this lower-bound on the specific Wilson coefficient will further reduce when more than one coupling are triggered simultaneously. However, the rest of our analysis is focused on the scenario where the sole operator contributes to the relic density.
For numerical computation of the DM relic density, we have used MadDM [81,82], which implements the exact and closed expression for thermal averaged annihilation crosssection as discussed by the authors in the reference [79]. The input model file required by the MadDM is generated using FeynRules [83,84], which calculates all the required couplings and Feynman rules by using the full Lagrangian given in equations (3), (6) and (8) corresponding to Majorana, scalar and vector DM particles respectively. We have further verified and validated our numerical results from MicrOMEGAs [85]. The constant relic density contours are depicted in figures 1, 2, and 3 corresponding to Majorana, real scalar and real vector DM candidates respectively in the plane defined by The thermal averaged DM pair annihilation cross-section is found to be ∼ 2 × 10 −26 cm 3 s −1 corresponding to a relic density of 0.1198, which essentially provides the functional dependence and shape-profile of the Wilson coefficients with the varying DM mass corresponding to each four-point effective operator. In order to understand the shape profile of the relic density contours, σ v may be expanded analytically in power-series of |v DM | 2 and then the contribution of the leading terms may be investigated. The velocity of the DM at freeze-out is considered to be 0.3 c. We observe that • due to the constrained thermal averaged cross-section, the cut-off Λ for Majorana DM scalar and pseudo-scalar operators in figures 1a and 1b varies as (m f m χ ) 1/3 and therefore shows a steep rise for top-philic DM in comparison to the bottom-philic case. However, the pseudo-scalar contribution dominates over its scalar counterpart due to its contribution from the leading velocity independent term in equation (B1b).
The cut-offs for leading contributions from the quark-philic twist operator and velocity suppressed gluon-philic scalar, pseudo-scalar, and twist operators are found to increase monotonously with increasing DM mass as m 3/4 χ .
• the s-wave contribution in the thermal averaged annihilation cross-section σ AV | v χ | given in equation (B1c) is proportional to m 2 f which follows from the chirality conserving property of axial-vector operator [87]. The variation of σ AV | v χ | with the increasing DM mass is found to be largely determined by the terms proportional to m 2 χ in the converging power-series expansion as explicitly shown in equation (B1c). As a consequence, to satisfy the relic density constraint, the cut-off varies with respect to DM mass as (m χ | − → v χ |) 1/2 , i.e, faster than scalar and pseudo-scalar but slower than the twist operator, as shown in figures 1a and 1b.
• the cut-offs for the scalar and twist operators induced by heavy-quark-philic scalar and vector DM candidates are depicted in figures 2a, 2b and 3a, 3b respectively. In case of scalar interactions, the cuts-offs are found to be independent of the DM masses, whereas the cut-off for the heavy-quark-philic vector DM pseudo-scalar operator is found to be monotonously increasing as (m V 0 ) 1/2 . In case of twist operators, the cut-off variation goes as (m φ 0 ) 1/2 for scalar DM, while for vector DM scenario, it is d-wave suppressed and goes as • the gluon-philic scalar, pseudo-scalar, and twist interactions corresponding to Majorana DM are p-wave suppressed and the cut-off dependences are depicted in figure 1c.
The cut-offs corresponding to gluon-philic scalar, pseudo-scalar, and twist operators for vector DM candidates are constrained to vary as ( respectively, as shown in figure 3c. The observed relative suppression of the constrained cut-off for the twist operator is due to the velocity dependence.
The same suit follows for the scalar DM candidates corresponding to the gluon-philic scalar and twist operators in figure 2c.

B. Indirect detection of DM pairs
Today, the WIMP Dark Matter in the universe is expected to be trapped in large gravitational potential wells, which further enhances the number density of DM in the region, resulting in frequent collisions among themselves. This facilitates DM-DM pair annihilation into a pair of SM particles (photons, leptons, hadronic jets, etc) at the Galactic centre, in Dwarf Spheroidal galaxies (dSphs), Galaxy clusters, and Galactic halos. The dwarf spheroidal satellite galaxies of the Milky Way are especially promising targets for DM indirect detection due to their large dark matter content, low diffuse Galactic γ-ray foregrounds as they travel the galactic distance, and lack of conventional astrophysical γray production mechanisms. Their flux is observed by the satellite-based γ-ray observatory  produced by DM pair annihilation hadronize partially to neutral pions, which then decay to photon pairs (π 0 → γγ) with a 99% branching fraction, yielding a broad spectrum of photons. In order to compare the photon spectra observed in FermiLAT [15] and H.E.S.S [17] experiments, the photon flux resulting from the DM pair annihilations is realised by inter-facing the MadDM algorithm [82] with the showering and hadronization simulated using PYTHIA 8.0 code [96]. In this subsection, we compare the thermally averaged DM annihilation cross-sections corresponding to the bb, tt and gg channels in the dSphs environment with the upper bounds calculated from the respective recasted experimental limits of the observed photon spectra.
The analytic expressions of thermally averaged DM pair annihilation cross-sections with what is known for other fermions in the literature [40,46,97,98]. These are derived from the cross-sections in the Appendix A in which the center of mass energy squared s It is to be noted that DM velocity | v DM | is roughly of the order of ∼ 10 −3 c at the center of Galaxy and 10 −5 c at dSphs in contrast to that of 10 −1 c at freeze-out. In comparison to all chiral blind s-dominated processes, the leading term contributions from the p-wave and d-wave channels We observe that the shape profile of the axial-vector induced contribution to the thermal averaged DM pair annihilation cross-section with varying DM mass is governed by the   The Majorana, real scalar, and real vector DM-nucleon scattering cross-sections driven by the heavy-quark scalar, axial-vector, and twist-2 currents are given as is the reduced mass for the respective DMnucleon system and the scale µ is taken to be Z 0 -boson mass. For the computation of scale dependent α s (µ) , α s (Λ) and g (x; Λ), we access CTEQ6l1 [99] PDF data set from LHAPDF6 [100] library. However, the scale dependence of the Wilson coefficients is found to be smaller than that of α s , as noted by the authors of the reference [101]. It is to be noted that the observed logarithmically enhanced one-loop induced scattering cross-sections in equations (12c), (13b) and (14c) result from the explicit momentum dependence in the twist-2 operator interaction Lagrangian given in equations (3), (6) and (8) The analytical expressions for the gluon-philic DM-nucleon scattering cross-sections are in agreement with those given in reference [74].
We display the spin independent Majorana, real scalar and real vector DM-nucleon scat- extracted by interpolating the data for event rates from the experimental data in references [102] and [103] for m DM less than 1 TeV.

Nuclear Recoil Spectrum
For a fixed target detector exposure T and target nucleus mass m nuc. , the differential nuclear recoil event rate w.r.t. nuclear recoil energy dR nuc. /dE r corresponding to the DM-nucleon scattering cross-section σ N and DM velocity | v| is given as where F (| q| r n ) is the Helm form factor taking account of the non-vanishing finite size of the nucleus [104]. Here r n ≡ 1.2 × A 1/3 is the effective radius of the nucleus with atomic mass A and | q| is the momentum transfer corresponding to the recoil energy E r . Following reference [105], the velocity integral for normalized Maxwellian DM velocity distribution is solved as For numerical computation we have taken Earth's velocity relative to galactic frame to be | v E | = 232 km/s, | v 0 | = 220 km/s and the escape velocity | v esc | = 544 km/s. Further, in order to incorporate the detector based effects, the nucleus recoil event rate dR nuc. / dE r is convoluted with the detector efficiency [6]. Integrating over the recoil energy from E th ∼ 4.9 KeV to the maximum E max for a fixed duration and size of the detector, we can estimate the expected recoil nucleus events.
As an illustration, we predict and plot the probable number of Xe nuclear recoil events

IV. SUMMARY AND CONCLUSIONS
In this article we assess the viability of the b-quark-philic, t-quark-philic, and gluon-philic self-conjugated spin 1/2, 0 and 1 DM candidates in the EFT approach. We have formulated the generalised effective interaction Lagrangian induced by the scalar, pseudo-scalar, axialvector, and twist-2 operators for the real particles in section II.
For a given DM mass, the relic abundance of Majorana, scalar, and vector DM is computed in section III A using the thermally averaged cross-sections in appendix B. Figures 1, 2 and 3 show twenty-eight interaction strengths in the form of Wilson coefficients C q, g DM i O j / Λ n in TeV −n (eleven for Majorana DM, six for scalar DM and eleven for vector DM) satisfying the relic density constraint Ω DM h 2 ≈ 0.1198 [4].
Using the constrained Wilson coefficients, we study the thermal averaged DM pair annihilation cross-sections for the indirect detection of varying DM masses (0.01 -2 TeV) in figures 4, 5 and 6. The contributions driven by the lower-limit of the effective couplings to the annihilation cross-sections in bb, tt and g g channels are found to be consistent when compared with the upper limits of the bb annihilation cross-sections obtained from FermiLAT [15] and H.E.S.S. [17] indirect experiments.
The scattering of the incident heavy-quark-philic DM particle off the static nucleon in- coefficients corresponding to scalar and twist-2 type-2 operators from the relic density constraint [4] and direct detection limits from XENON-1T experiment [6], shown in red, blue and green for the b-quark, t-quark and gluon-philic scalar DM respectively.
As mentioned in the text, we include the study of gluon-philic Majorana DM interactions induced by the scalar, pseudo-scalar and twist-2 type-1 operators. The annihilation crosssections of the Majorana DM pair to a pair of gluons are given as The annihilation cross-sections of a pair of real spin 0 DM φ 0 to a pair of third generation heavy quarks induced by the scalar and twist-2 type-2 operators are given as: The cross-sections for the pair of scalar DM annihilation into a pair of gluons induced by scalar and twist-2 type-2 operators are given as Similarly, the production of a pair of third generation heavy quarks as a result of the annihilation of a pair of real vectors DM induced by the scalar, pseudo-scalar, axial-vector, and twist-2 type-2 operators is given as The annihilation cross-sections of a pair of vector DM induced by the scalar and twist-2 type-2 gluon currents are given as The thermal average of the heavy-quark-philic Majorana DM pair annihilation crosssections given in equations (A1a), (A1b), (A1c) and (A1d) corresponding to scalar and pseudo-scalar operators respectively are given as The axial-vector contribution is analytically expanded in terms of two independently converging series expansions in terms of | v χ | 2 and given as And finally the contribution from the twist-2 type-1 operator induced by the heavy quarkphilic Majorana DM interaction is given as Using annihilation cross-sections in (A2a), (A2b) and (A2c) the thermal averaged annihilation cross-sections for the gluon-philic Majorana DM are given as The thermal average of heavy-quark-philic scalar DM pair annihilation cross-sections dipalyed in equations (A3a) and (A3b) are given as Similarly, the thermal average of heavy-quark-philic scalar DM pair annihilation crosssections displayed in equations (A4a) and (A4b) are given as The thermal average of the vector DM pair annihilation cross-sections given in equations (A5a), (A5b), (A5c) and (A5d) corresponding to the scalar, pseudo-scalar, axial-vector and twist-2 type-2 operators, respectively, are given as Similarly, the thermal average of the vector DM pair annihilation to gluon channels as displayed in equations (A6a), (A6b) and (A6c) corresponding to the scalar, pseudo-scalar and twist-2 currents respectively are given as The dimensionless one-loop integrals I gg S , I gg AV , and I gg T are defined as where the square of the four momentum transferred q 2 ≈ 2 m N E r (E r 100 KeV is the recoil energy of the nucleon). m Q and m N are the concerned heavy quark mass (m b / m t ) running in the loop and nucleon mass respectively.
We observe that the one-loop effective contact interactions generated from quark scalar, twist and axial vector currents contain the scalar, pseudo-scalar and twist gluon operators respectively. We perform the non-relativistic reduction of these operators for zero momentum partonic gluons and evaluate the gluonic operators between the nucleon states. The zero momentum gluon contribution to the hadronic matrix element f N TG ≈ .923 is extracted in terms of light quark as [106][107][108] According to [108] and [109], the pseudoscalar gluon operator between nucleon states is s , while the contribution from heavy quarks is found to be vanishingly small [110]. It is important to observe that a quark axial-vector current is related to the gluonic pseudoscalar current by PCAC [109]. The zero momentum nucleonic matrix element for gluon twist-2 operators is defined as where g (2; µ R ) =