Nuclear structure factors for general spin-independent WIMP-nucleus scattering

TL;DR

This paper provides detailed nuclear structure factors for spin-independent WIMP-nucleus scattering, incorporating advanced nuclear models and two-nucleon currents, to improve dark matter detection predictions.

Contribution

It introduces comprehensive structure factors including spin-2 operators and two-nucleon currents using large-scale shell model calculations, expanding previous work on xenon to other targets.

Findings

Results for fluorine, silicon, argon, and germanium nuclei.

Inclusion of spin-2 operator contributions.

Public availability of results in a Python notebook.

Abstract

We present nuclear structure factors that describe the generalized spin-independent coupling of weakly interacting massive particles (WIMPs) to nuclei. Our results are based on state-of-the-art nuclear structure calculations using the large-scale nuclear shell model. Starting from quark- and gluon-level operators, we consider all possible coherently enhanced couplings of spin-1/2 and spin-0 WIMPs to one and two nucleons up to third order in chiral effective field theory. This includes a comprehensive discussion of the structure factors corresponding to the leading two-nucleon currents covering, for the first time, the contribution of spin-2 operators. We provide results for the most relevant nuclear targets considered in present and planned dark matter direct detection experiments: fluorine, silicon, argon, and germanium, complementing our previous work on xenon. All results are also publicly available in a Python notebook.

Figures

Tables

Table 1. Table 1: Nuclear structure factors and associated QCD and NR operators. Isospin indices are suppressed because the identification pertains to isoscalar and isovector combinations alike. The table summarizes all operators that can be coherently enhanced. Likewise, the table only shows coherently enhanced contributions for operators with velocity v ⟂ superscript 𝑣 perpendicular-to v^{\perp} , which in nuclei produce both a term that behaves as q / m N 𝑞 subscript 𝑚 𝑁 q/m_{N} and a remainder determined by the target velocity v T ⟂ ∼ 10 − 3 similar-to subscript superscript 𝑣 perpendicular-to 𝑇 superscript 10 3 v^{\perp}_{T}\sim 10^{-3} (see Appendix  A for more details). For comparison the last rows for S i ​ j subscript 𝑆 𝑖 𝑗 S_{ij} show operators that lead to SD interactions. In the chiral counting we have M π = 𝒪 ​ ( p ) subscript 𝑀 𝜋 𝒪 𝑝 M_{\pi}=\mathcal{O}(p) , v ⟂ = 𝒪 ​ ( p 2 ) superscript 𝑣 perpendicular-to 𝒪 superscript 𝑝 2 v^{\perp}=\mathcal{O}(p^{2}) , with relativistic corrections counted as ∂ / m N = 𝒪 ​ ( p 2 ) subscript 𝑚 𝑁 𝒪 superscript 𝑝 2 \partial/m_{N}=\mathcal{O}(p^{2}) and ∂ / m χ = 𝒪 ​ ( p 2 ) subscript 𝑚 𝜒 𝒪 superscript 𝑝 2 \partial/m_{\chi}=\mathcal{O}(p^{2}) . Two-body structure factors, ℱ π subscript ℱ 𝜋 \mathcal{F}_{\pi} and ℱ b subscript ℱ b \mathcal{F}_{\text{b}} , cannot be matched onto single-nucleon NREFT operators, and the quasi-coherence of ℱ Φ ′′ superscript ℱ superscript Φ ′′ \mathcal{F}^{\Phi^{\prime\prime}} is characterized by ξ ∼ 0.2 similar-to 𝜉 0.2 \xi\sim 0.2 . Finally, the last column indicates whether each entry interferes with the leading 𝒪 1 subscript 𝒪 1 \mathcal{O}_{1} operator.

Structure factor	QCD operators	Chiral scaling	NR operators	Overall scaling	Interference with ${𝒪}_1$
${ℱ}^M$	$\overline{\chi }\chi m_q\overline{q}q$	$𝒪\left(M_{\pi }^2\right)=𝒪\left(p^2\right)$	${𝒪}_1$	$𝒪\left(\frac{M_{\pi }^2}{{\Lambda }_{\chi }^2}A\right)$	yes
	$\overline{\chi }\chi {\theta }_{\mu }^{\mu }$	$𝒪\left(1\right)$	${𝒪}_1$	$𝒪\left(A\right)$	yes
	$\overline{\chi }{\gamma }_{\mu }i{\partial }_{\nu }\chi {\overline{\theta }}^{\mu \nu }$	$𝒪\left(1\right)$	${𝒪}_1$	$𝒪\left(A\right)$	yes
	$\overline{\chi }{\gamma }^{\mu }\chi \overline{q}{\gamma }_{\mu }q$	$𝒪\left(1\right)$	${𝒪}_1$	$𝒪\left(A\right)$	yes
		$𝒪\left(\frac{q{v}^{⟂}}{m_{\chi }+m_N}\right)=𝒪\left(p^4\right)$	${𝒪}_5$	$𝒪\left(\frac{q{v}_T^{⟂}}{m_{\chi }+m_N}A\right)$	no
	$\overline{\chi }{\sigma }^{\mu \nu }\chi \overline{q}{\sigma }_{\mu \nu }q$	$𝒪\left(\frac{q^2}{m_N{m}_{\chi }}\right)=𝒪\left(p^4\right)$	${𝒪}_1$	$𝒪\left(\frac{q^2}{m_N{m}_{\chi }}A\right)$	yes
		$𝒪\left(\frac{q{v}^{⟂}}{m_N}\right)=𝒪\left(p^4\right)$	${𝒪}_5$	$𝒪\left(\frac{q{v}_T^{⟂}}{m_N}A\right)$	no
	$\overline{\chi }{\gamma }^{\mu }{\gamma }_5\chi \overline{q}{\gamma }_{\mu }q$	$𝒪\left(v^{⟂}\right)=𝒪\left(p^2\right)$	${𝒪}_8$	$𝒪\left(v_T^{⟂}A\right)$	no
	$\overline{\chi }{\sigma }^{\mu \nu }i{\gamma }_5\chi \overline{q}{\sigma }_{\mu \nu }q$	$𝒪\left(\frac{q}{m_N}\right)=𝒪\left(p^2\right)$	${𝒪}_{11}$	$𝒪\left(\frac{q}{m_N}A\right)$	no
	$\overline{\chi }i{\gamma }_5\chi m_q\overline{q}q$	$𝒪\left(\frac{q{M}_{\pi }^2}{m_{\chi }}\right)=𝒪\left(p^4\right)$	${𝒪}_{11}$	$𝒪\left(\frac{q{M}_{\pi }^2}{m_{\chi }{\Lambda }_{\chi }^2}A\right)$	no
${ℱ}_{\pi }$	$\overline{\chi }\chi m_q\overline{q}q$	$𝒪\left(M_{\pi }^3\right)=𝒪\left(p^3\right)$	–	$𝒪\left(\frac{M_{\pi }^3}{{\Lambda }_{\chi }^3}A\right)$	yes
	$\overline{\chi }\chi {\theta }_{\mu }^{\mu }$	$𝒪\left(M_{\pi }^3\right)=𝒪\left(p^3\right)$	–	$𝒪\left(\frac{M_{\pi }^3}{{\Lambda }_{\chi }^3}A\right)$	yes
	$\overline{\chi }{\gamma }_{\mu }i{\partial }_{\nu }\chi {\overline{\theta }}^{\mu \nu }$	$𝒪\left(M_{\pi }^3\right)=𝒪\left(p^3\right)$	–	$𝒪\left(\frac{M_{\pi }^3}{{\Lambda }_{\chi }^3}A\right)$	yes
${ℱ}_{\text{b}}$	$\overline{\chi }\chi {\theta }_{\mu }^{\mu }$	$𝒪\left(M_{\pi }^3\right)=𝒪\left(p^3\right)$	–	$𝒪\left(\frac{M_{\pi }^3}{{\Lambda }_{\chi }^3}A\right)$	yes
	$\overline{\chi }{\gamma }_{\mu }i{\partial }_{\nu }\chi {\overline{\theta }}^{\mu \nu }$	$𝒪\left(M_{\pi }^3\right)=𝒪\left(p^3\right)$	–	$𝒪\left(\frac{M_{\pi }^3}{{\Lambda }_{\chi }^3}A\right)$	yes
${ℱ}^{{\Phi }^{\prime \prime }}$	$\overline{\chi }{\gamma }^{\mu }\chi \overline{q}{\gamma }_{\mu }q$	$𝒪\left(q{v}^{⟂}\right)=𝒪\left(p^3\right)$	${𝒪}_3$	$𝒪\left(\frac{q^2}{m_N{\Lambda }_{\chi }}\xi A\right)$	yes
	$\overline{\chi }{\sigma }^{\mu \nu }\chi \overline{q}{\sigma }_{\mu \nu }q$	$𝒪\left(\frac{q{v}^{⟂}}{m_{\chi }}\right)=𝒪\left(p^4\right)$	${𝒪}_3$	$𝒪\left(\frac{q^2}{m_N{m}_{\chi }}\xi A\right)$	yes
$S_{ij}$	$\overline{\chi }{\gamma }^{\mu }{\gamma }_5\chi \overline{q}{\gamma }_{\mu }{\gamma }_{5}q$	$𝒪\left(1\right)$	${𝒪}_4$, ${𝒪}_6$	$𝒪\left(1\right)$	no
	$\overline{\chi }{\sigma }^{\mu \nu }\chi \overline{q}{\sigma }_{\mu \nu }q$	$𝒪\left(1\right)$	${𝒪}_4$	$𝒪\left(1\right)$	no
	$\overline{\chi }{\gamma }_5\chi \overline{q}{\gamma }_{5}q$	$𝒪\left(1\right)$	${𝒪}_6$	$𝒪\left(1\right)$	no

Table 2. Table 2: Root-mean-square charge radii ( ⟨ r 2 ⟩ ch 1 / 2 superscript subscript delimited-⟨⟩ superscript 𝑟 2 ch 1 2 \langle r^{2}\rangle_{\text{ch}}^{1/2} ) of ground states (gs), electric quadrupole ( Q 𝑄 Q ) and magnetic dipole ( μ 𝜇 \mu , in units of nuclear magnetons, n.m.) moments of ground and lowest-excited states, and nuclear matrix elements of selected electric quadrupole [B(E2)] and magnetic dipole [B(M1)] transitions between low-lying states of stable fluorine and silicon isotopes. Theoretical charge radii use calculated orbital occupancies combined with the harmonic oscillator length b 𝑏 b from Eq. ( 64 ). Electromagnetic moments and transitions are obtained with effective neutron and proton electric charges e n = 0.45 subscript 𝑒 𝑛 0.45 e_{n}=0.45 and e p = 1.36 subscript 𝑒 𝑝 1.36 e_{p}=1.36 , and effective orbital ( l 𝑙 l ) and spin ( s 𝑠 s ) g 𝑔 g -factors g p l = 1.16 superscript subscript 𝑔 𝑝 𝑙 1.16 g_{p}^{l}=1.16 , g p s = 5.15 superscript subscript 𝑔 𝑝 𝑠 5.15 g_{p}^{s}=5.15 for protons and g n l = − 0.09 superscript subscript 𝑔 𝑛 𝑙 0.09 g_{n}^{l}=-0.09 , g n s = − 3.55 superscript subscript 𝑔 𝑛 𝑠 3.55 g_{n}^{s}=-3.55 for neutrons  usdb_obs . Quadrupole moments Q > 0 𝑄 0 Q>0 ( Q < 0 𝑄 0 Q<0 ) indicate nuclei with prolate (oblate) deformation. Theoretical results are compared to experimental data from Refs.  Angeli:2013epw ; nndc ; Raghavan:1989zz .

Nucleus	State /	${⟨r^2⟩}_{\text{ch}}^{1/2}$ [fm]		$Q$ [$e$fm2]		$\mu$ [n.m.]		B(E2) [$e^2$fm4]		B(M1) [n.m.2]
	Transition	Th	Exp	Th	Exp	Th	Exp	Th	Exp	Th	Exp
19F	$1/2_{\text{gs}}^{+}$	$2.83$	$2.898\left(2\right)$	–	–	$2.9$	$2.628868\left(8\right)$
	$5/2_1^{+}\to 1/2_{\text{gs}}^{+}$							$19$	$20.9\left(2\right)$	–	–
	$3/2_1^{+}\to 5/2_1^{+}$									$3.4$	$4.1\left(2.5\right)$
	$9/2_1^{+}\to 5/2_1^{+}$							$19$	$25\left(3\right)$	–	–
28Si	$0_{\text{gs}}^{+}$	$3.19$	$3.122\left(2\right)$	–	–	–	–
	$2_1^{+}$			$+19$	$+16\left(3\right)$	$1.1$	$1.12\left(18\right)$
	$2_1^{+}\to 0_{\text{gs}}^{+}$							$67$	$70\left(3\right)$	–	–
	$4_1^{+}\to 2_1^{+}$							$110$	$87\left(10\right)$	–	–
	$0_2^{+}\to 2_1^{+}$							$82$	$50\left(3\right)$	–	–
	$3_1^{+}\to 2_1^{+}$							$6\times {10}^{-4}$	$0.007\left(2\right)$	$2\times {10}^{-4}$	$48\left(4\right)\times {10}^{-5}$
	$0_3^{+}\to 2_1^{+}$							$0.6$	$1.4\left(1\right)$	–	–
29Si	$1/2_{\text{gs}}^{+}$	$3.20$	$3.118\left(5\right)$	–	–	$-0.49$	$-0.55529\left(3\right)$
	$3/2_1^{+}\to 1/2_{\text{gs}}^{+}$							$31$	$22\left(2\right)$	$0.011$	$0.063\left(2\right)$
	$3/2_2^{+}\to 1/2_{\text{gs}}^{+}$							$33$	$29\left(12\right)$	$0.21$	$0.116\left(7\right)$
	$5/2_2^{+}\to 3/2_1^{+}$							$64$	$47\left(7\right)$	$0.20$	$0.15\left(1\right)$
	$5/2_2^{+}\to 5/2_1^{+}$							$2.5$	$4\left(4\right)$	$0.18$	$0.18\left(3\right)$
	$7/2_1^{+}\to 3/2_1^{+}$							$50$	$40\left(12\right)$	–	–
	$7/2_1^{+}\to 5/2_1^{+}$							$1.8$	$1.0\left(6\right)$	$0.013$	$0.04\left(1\right)$
30Si	$0_{\text{gs}}^{+}$	$3.21$	$3.134\left(4\right)$	–	–	–	–
	$2_1^{+}$			$+2.1$	$-5\left(6\right)$	$0.70$	$0.76\left(18\right)$
	$2_1^{+}\to 0_{\text{gs}}^{+}$							$48$	$47\left(6\right)$	–	–
	$2_2^{+}\to 0_{\text{gs}}^{+}$							$14$	$33\left(9\right)$	–	–
	$2_2^{+}\to 2_1^{+}$							$81$	$50\left(30\right)$	$0.18$	$0.16\left(5\right)$
	$1_1^{+}\to 0_{\text{gs}}^{+}$							–	–	$0.007$	$0.009\left(3\right)$
	$1_1^{+}\to 2_1^{+}$							$2.3$	$8\left(6\right)$	$0.29$	$0.16\left(4\right)$

Table 3. Table 3: Same as Table  2 , for stable argon and germanium isotopes. Here electromagnetic moments and transitions are obtained with effective neutron and proton electric charges e n = 0.5 subscript 𝑒 𝑛 0.5 e_{n}=0.5 , and e p = 1.5 subscript 𝑒 𝑝 1.5 e_{p}=1.5 , and bare g 𝑔 g -factors  Caurier:2004gf . Quadrupole moments Q > 0 𝑄 0 Q>0 ( Q < 0 𝑄 0 Q<0 ) indicate nuclei with prolate (oblate) deformation. Theoretical results are compared to experimental data from Refs.  Angeli:2013epw ; nndc ; Raghavan:1989zz .

Nucleus	State /	${⟨r^2⟩}_{\text{ch}}^{1/2}$ [fm]		$Q$ [$e$fm2]		$\mu$ [n.m.]		B(E2) [$e^2$fm4]		B(M1) [n.m.2]
	Transition	Th	Exp	Th	Exp	Th	Exp	Th	Exp	Th	Exp
40Ar	$0_{\text{gs}}^{+}$	$3.43$	$3.427\left(3\right)$	–	–	–	–
	$2_1^{+}$			$+2.6$	$+1\left(4\right)$	$-0.54$	$-0.04\left(6\right)$
	$2_1^{+}\to 0_{\text{gs}}^{+}$							$50$	$73\left(3\right)$	–	–
	$0_2^{+}\to 2_1^{+}$							$29$	$43\left(7\right)$	–	–
	$2_2^{+}\to 0_{\text{gs}}^{+}$							$0.7$	$10\left(2\right)$	–	–
	$2_2^{+}\to 2_1^{+}$							$55$	$150\left(50\right)$	$0.016$	$0.07\left(1\right)$
	$4_1^{+}\to 2_1^{+}$							$36$	$43\left(8\right)$	–	–
	$6_1^{+}\to 4_1^{+}$							$16$	$13.6\left(5\right)$	–	–
70Ge	$0_{\text{gs}}^{+}$	$4.05$	$4.0414\left(12\right)$	–	–	–	–
	$2_1^{+}$			$+23$	$+4\left(3\right)$	$0.96$	$0.91\left(5\right)$
	$2_1^{+}\to 0_{\text{gs}}^{+}$							$240$	$360\left(7\right)$	–	–
	$0_2^{+}\to 2_1^{+}$							$36$	$820\left(120\right)$	–	–
	$2_2^{+}\to 0_{\text{gs}}^{+}$							$8.0$	$9\left(1\right)$	–	–
	$2_2^{+}\to 2_1^{+}$							$16$	$1100\left(190\right)$	$0.022$	$0.003\left(2\right)$
	$2_2^{+}\to 0_2^{+}$							$270$	$270\left(50\right)$	–	–
	$4_1^{+}\to 2_1^{+}$							$370$	$430\left(90\right)$	–	–
72Ge	$0_{\text{gs}}^{+}$	$4.07$	$4.0576\left(12\right)$	–	–	–	–
	$2_1^{+}$			$+16$	$-13\left(6\right)$	$0.55$	$0.77\left(5\right)$
	$2_1^{+}\to 0_{\text{gs}}^{+}$							$260$	$418\left(7\right)$	–	–
	$2_1^{+}\to 0_2^{+}$							$60$	$317\left(5\right)$	–	–
	$2_2^{+}\to 0_{\text{gs}}^{+}$							$29$	$2.3\left(4\right)$	–	–
	$2_2^{+}\to 0_2^{+}$							$15$	$0.5\left(1\right)$	–	–
	$2_2^{+}\to 2_1^{+}$							$360$	$1100\left(180\right)$	$0.023$	$29\left(9\right)\times {10}^{-5}$
	$4_1^{+}\to 2_1^{+}$							$430$	$660\left(90\right)$	–	–
73Ge	$9/2_{\text{gs}}^{+}$	$4.07$	$4.0632\left(14\right)$	$-15$	$-17\left(3\right)$	$-1.0$	$-0.8794677\left(2\right)$
	$5/2_1^{+}$			$+20$	$+70\left(8\right)$	$-0.90$	$-1.08\left(3\right)$
	$7/2_1^{+}\to 9/2_{\text{gs}}^{+}$							$260$	$114\left(7\right)$	$0.002$	$14\left(9\right)\times {10}^{-4}$
	$7/2_2^{+}\to 9/2_{\text{gs}}^{+}$							$26$	$740\left(140\right)$	$0.014$	$57\left(4\right)\times {10}^{-4}$
	$13/2_1^{+}\to 9/2_{\text{gs}}^{+}$							$250$	$540\left(40\right)$	–	–
74Ge	$0_{\text{gs}}^{+}$	$4.08$	$4.0742\left(12\right)$	–	–	–	–
	$2_1^{+}$			$+0.7$	$-19\left(2\right)$	$0.49$	$0.87\left(4\right)$
	$2_1^{+}\to 0_{\text{gs}}^{+}$							$310$	$609\left(7\right)$	–	–
	$2_2^{+}\to 0_{\text{gs}}^{+}$							$3.2$	$13\left(2\right)$	–	–
	$2_2^{+}\to 2_1^{+}$							$470$	$790\left(110\right)$	$4\times {10}^{-5}$	$18\left(3\right)\times {10}^{-4}$
	$4_1^{+}\to 2_1^{+}$							$450$	$760\left(60\right)$	–	–
76Ge	$0_{\text{gs}}^{+}$	$4.09$	$4.0811\left(12\right)$	–	–	–	–
	$2_1^{+}$			$-14$	$-19\left(6\right)$	$0.42$	$0.84\left(5\right)$
	$2_1^{+}\to 0_{\text{gs}}^{+}$							$300$	$550\left(19\right)$	–	–
	$2_2^{+}\to 0_{\text{gs}}^{+}$							$1.2$	$17\left(4\right)$	–	–
	$2_2^{+}\to 2_1^{+}$							$400$	$800\left(170\right)$	$0.006$	$14\left(7\right)\times {10}^{-4}$
	$4_1^{+}\to 2_1^{+}$							$430$	$730\left(170\right)$	–	–

Table 4. Table 4: Theoretical root-mean-square charge radii of stable xenon isotopes compared to experimental data from Ref.  Angeli:2013epw .

		128Xe	129Xe	130Xe	131Xe
${⟨r^2⟩}_{\text{ch}}^{1/2}$	Th	$4.75$	$4.75$	$4.76$	$4.77$
$\left[\text{fm}\right]$	Exp	$4.7774\left(50\right)$	$4.7775\left(50\right)$	$4.7818\left(49\right)$	$4.7808\left(49\right)$
		132Xe	134Xe	136Xe
${⟨r^2⟩}_{\text{ch}}^{1/2}$	Th	$4.77$	$4.78$	$4.79$
$\left[\text{fm}\right]$	Exp	$4.7808\left(49\right)$	$4.7899\left(47\right)$	$4.7964\left(47\right)$

Table 5. Table 5: C S subscript 𝐶 𝑆 C_{S} values, renormalized to the nuclear binding energy, for all isotopes considered in this work.

	19F	28Si	29Si	30Si	40Ar	70Ge
$C_S\left[{\text{GeV}}^{-2}\right]$	$-68.9$	$-75.2$	$-75.0$	$-75.0$	$-76.6$	$-85.6$
	72Ge	73Ge	74Ge	76Ge	128Xe	129Xe
$C_S\left[{\text{GeV}}^{-2}\right]$	$-85.9$	$-85.9$	$-86.2$	$-86.4$	$-98.0$	$-98.2$
	130Xe	131Xe	132Xe	134Xe	136Xe
$C_S\left[{\text{GeV}}^{-2}\right]$	$-98.5$	$-98.6$	$-98.8$	$-99.1$	$-99.3$

Table 6. Table 6: Maximal power m 𝑚 m in Eq. ( 63 ) for each structure factor and nuclear target considered in this work.

Structure factor	F	Si	Ar	Ge	Xe
${ℱ}_{\pm }^M$	2	2	3	4	5
${ℱ}_{\pm }^{{\Phi }^{\prime \prime }}$	1	1	2	3	4
${ℱ}_{\pi }$	3	3	3	4	5
${ℱ}_{\text{b}}$	2	2	3	4	5

Equations

$$\frac{\text{d}\sigma }{\text{d}{q}^2}$$

$$\displaystyle+c_{\text{b}}\mathcal{F}_{\text{b}}(q^{2})+\frac{q^{2}}{2m_{N}^{2}}\sum_{I=\pm}c_{I}^{\Phi^{\prime\prime}}\mathcal{F}_{I}^{\Phi^{\prime\prime}}(q^{2})\bigg{|}^{2}$$

$$\displaystyle+\frac{1}{4\pi v^{2}}\sum_{i=5,8,11}\bigg{|}\sum_{I=\pm}\xi_{i}(q,v^{\perp}_{T})c_{I}^{M,i}\mathcal{F}_{I}^{M}(q^{2})\bigg{|}^{2},$$

$${\xi }_5=\frac{{\mu }_{N}q{v}_T^{\perp }}{2m_{\chi }{m}_N},{\xi }_8=v_T^{\perp },{\xi }_{11}=-\frac{q}{2m_{\chi }},$$

$${\mathcal{L}}_{\chi }$$

$${\mathcal{L}}_{\chi }^{(6)}$$

$$+C_q^{AV}\bar{\chi }{\gamma }^{\mu }{\gamma }_5\chi \bar{q}{\gamma }_{\mu }q$$

$$\displaystyle+C_{q}^{TT}\bar{\chi}\sigma^{\mu\nu}\chi\,\bar{q}\sigma_{\mu\nu}q+\tilde{C}_{q}^{TT}\bar{\chi}\sigma^{\mu\nu}i\gamma_{5}\chi\,\bar{q}\sigma_{\mu\nu}q\Big{]},$$

$${\mathcal{L}}_{\chi }^{(7)}$$

$$\displaystyle+\sum_{q}\Big{(}C_{q}^{PS}+\frac{8\pi}{9}\tilde{C}^{\prime S}_{g}\Big{)}\bar{\chi}i\gamma_{5}\chi\,m_{q}\bar{q}q$$

$$\displaystyle-\frac{8\pi}{9}C^{\prime S}_{g}\bar{\chi}\chi\,\theta^{\mu}_{\mu}-\frac{8\pi}{9}\tilde{C}^{\prime S}_{g}\bar{\chi}i\gamma_{5}\chi\,\theta^{\mu}_{\mu}\bigg{\}},$$

$${\mathcal{L}}_{\chi }^{(8)}$$

$$C_g^{\prime S}$$

$${\tilde{C}}_g^{\prime S}$$

$${\theta }_{\mu }^{\mu }=\sum _q{m}_q\bar{q}q-\frac{9}{8\pi }{\alpha }_s{G}_{\mu \nu }^a{G}_a^{\mu \nu }+\mathcal{O}({\alpha }_s^2).$$

$${\bar{\theta }}_q^{\mu \nu }$$

$${\bar{\theta }}_g^{\mu \nu }$$

$$c_{\pm }^M$$

$${\dot{c}}_{\pm }^M$$

$$\displaystyle+\frac{1}{4m_{N}^{2}}\big{(}f_{2}^{V,p}\pm f_{2}^{V,n}\big{)}+\frac{1}{2m_{\chi}m_{N}}\big{(}f_{1}^{T,p}\pm f_{1}^{T,n}\big{)}$$

$$\displaystyle+\frac{1}{m_{\chi}m_{N}}\big{(}f_{2}^{T,p}\pm f_{2}^{T,n}\big{)}\bigg{]},$$

$$c_{\pi }$$

$$c_{\pm }^{{\Phi }^{\prime \prime }}$$

$$\displaystyle+\frac{m_{N}}{m_{\chi}}\Big{(}1+\frac{\mu_{N}}{m_{N}}\Big{)}\big{(}f^{T,p}_{1}\pm f^{T,n}_{1}\big{)}\bigg{]},$$

$$c_{\pm }^{M,5}$$

$$\displaystyle+\Big{(}1+\frac{2m_{\chi}}{m_{N}}\Big{)}\big{(}f^{T,p}_{1}\pm f^{T,n}_{1}\big{)}+4\frac{m_{\chi}}{\mu_{N}}\big{(}f^{T,p}_{2}\pm f^{T,n}_{2}\big{)}\bigg{]},$$

$$c_{\pm }^{M,8}$$

$$c_{\pm }^{M,11}$$

$$\displaystyle-4\frac{m_{\chi}}{m_{N}}\big{(}\tilde{f}^{T,p}_{2}\pm\tilde{f}^{T,n}_{2}\big{)}\bigg{]},$$

$$f_{\pi }$$

$$f_{\pi }^{\theta }$$

$$f_{\pi }^{(2)}$$

$$f_u^{\pi }$$

$$f_d^{\pi }$$

$$f^{(2)}_{q,\pi}=\int_{0}^{1}\text{d}x\,x\big{(}q_{\pi}(x)+\bar{q}_{\pi}(x)\big{)},$$

$$\sum _q{f}_{q,\pi }^{(2)}+f_{g,\pi }^{(2)}=1.$$

$$f_{u,\pi }^{(2)}=f_{d,\pi }^{(2)}=⟨x{⟩}_{\text{v}}^{\pi }+2⟨x{⟩}_{\text{s}}^{\pi },f_{s,\pi }^{(2)}=2⟨x{⟩}_{\text{s}}^{\pi },$$

$$f_{u,\pi }^{(2)}$$

$$f_{s,\pi }^{(2)}$$

$$f_{u,\pi }^{(2)}=f_{d,\pi }^{(2)}=0.26(5),f_{s,\pi }^{(2)}=0.01(5),$$

$$f_{u,\pi }^{(2)}$$

$$f_{g,\pi }^{(2)}$$

$${\mathcal{L}}_{\chi }^{(5)}=\frac{C_F}{\Lambda }\bar{\chi }{\sigma }^{\mu \nu }\chi F_{\mu \nu }+\frac{{\tilde{C}}_F}{\Lambda }\bar{\chi }{\sigma }^{\mu \nu }\chi {\tilde{F}}_{\mu \nu },$$

$${\mathcal{M}}_{1,\text{NR}}^F$$

Full text

††preprint: INT-PUB-18-059

I Introduction

Astrophysical observations have established that more than three quarters of the matter content of the universe are composed of dark matter. The nature of dark matter, however, remains elusive, and its very existence is one of the most compelling pieces of evidence for physics beyond the Standard Model of particle physics. A key step to unveil the composition of dark matter would be its direct detection in the laboratory Baudis:2016qwx . Such endeavor is led by international collaborations that use atomic nuclei as targets, in their aim to detect the nuclear recoils resulting from the scattering of dark matter particles Angloher:2015ewa ; Amole:2016pye ; Armengaud:2016cvl ; Akerib:2016vxi ; Aprile:2017iyp ; Agnese:2017jvy ; Amaudruz:2017ekt ; Cui:2017nnn ; Agnese:2017njq ; Agnes:2018ves ; XMASS:2018bid ; Aprile:2018dbl . In spite of the very high sensitivities achieved by reduced backgrounds combined with extended exposures, there has been no conclusive evidence to date for the direct detection of dark matter. Next generations of experiments plan to push the frontier for dark matter direct detection by several orders of magnitude Aprile:2015uzo ; Akerib:2015cja ; PandaXxT ; Aalbers:2016jon , until the background from coherent neutrino–nucleus scattering Akimov:2017ade becomes dominant.

Direct detection experiments are motivated by extensions of the Standard Model that propose dark matter candidates interacting with quarks and gluons, the Standard Model fields that ultimately form atomic nuclei. Prominent candidates are weakly interacting massive particles (WIMPs) Bertone . Because the WIMPs forming dark matter would be non-relativistic (NR), their scattering off atomic nuclei would transfer energies and momenta much smaller than the nuclear or nucleon masses. As a consequence nucleons and nuclei, instead of quarks and gluons, become the relevant degrees of freedom. The interpretation of the present experimental limits and future values of WIMP–nucleus cross sections, therefore, naturally depends on the nuclear physics aspects of the scattering. This information is encoded in the so-called nuclear structure factors Engel:1992bf . While for some simple cases a phenomenological prescription can be a good approximation to the structure factor Lewin:1995rx ; Vietze:2014vsa , in general a good description of the nucleus, obtained with a dedicated nuclear many-body calculation, is needed.

In the absence of experimental data on the spin Baudis:2013bba ; McCabe:2015eia , momentum-transfer Rogers:2016jrx ; Fieguth:2018vob , or isospin dependence, the character of the WIMP–nucleus interaction is unknown. Nevertheless standard analyses usually assume so-called spin-independent (SI) interactions, which receive the coherent contribution of all nucleons in the nucleus. However, additional interactions are possible, and could be dominant if the SI coherence is compensated by suppressed values of the corresponding WIMP–nucleon couplings. Some examples analyzed experimentally are the so-called spin-dependent (SD) Amole:2016pye ; Aprile:2016swn ; Fu:2016ega ; Akerib:2017kat ; Aprile:2019dbj ; Amole:2019fdf or momentum-transfer-dependent interactions Angloher:2016jsl .

In order to organize different possible WIMP–nucleus interactions, two alternative schemes have been proposed recently. On the one hand, a non-relativistic effective field theory (NREFT) Fan:2010gt ; Fitzpatrick:2012ix ; Anand:2013yka based on the lowest-order operators that can describe the coupling of a WIMP to a nucleon. The NREFT approach proposes a set of one-body operators, considered in recent experimental analyses Schneck:2015eqa ; Aprile:2017aas ; Xia:2018qgs ; Angloher:2018fcs . On the other hand, an organization based on chiral effective field theory (ChEFT) Epelbaum:2008ga ; Machleidt:2011zz ; Hammer:2012id , a low-energy effective theory of quantum chromodynamics (QCD) that preserves the QCD symmetries, in particular capturing the important role played by pions at low energies. The ChEFT approach can be mapped onto the single-nucleon couplings of NREFT Hoferichter:2015ipa , implying certain interdependencies for the latter. In addition, ChEFT predicts consistent couplings of WIMPs to two nucleons Prezeau:2003sv ; Cirigliano:2012pq ; Menendez:2012tm ; Klos:2013rwa ; Cirigliano:2013zta ; Hoferichter:2016nvd ; Hoferichter:2017olk ; Andreoli:2018etf , reflecting that nucleons are strongly interacting in nuclei. Such contributions occur, e.g., when the WIMP couples to a virtual pion exchanged between the two nucleons, an effect recently constrained for the first time by direct-detection experiments Aprile:2018cxk , for the case of a scalar WIMP–quark interaction. Similar couplings to two nucleons are important in electromagnetic and weak transitions of atomic nuclei Gazit:2008ma ; Menendez:2011qq ; Bacca:2014tla ; Pastore:2017uwc . Furthermore, ChEFT provides a power counting that suggests a hierarchy, guided by QCD, for the expected importance of the different NREFT operators Hoferichter:2015ipa . This hierarchy is only tentative, because the couplings describing the interaction of the WIMP to the Standard Model fields are not known. A related ChEFT approach limited to WIMP interactions with one nucleon has been proposed in Refs. Bishara:2016hek ; Bishara:2017pfq .

In this work we follow Ref. Hoferichter:2016nvd and combine the ChEFT framework with large-scale shell model nuclear many-body calculations to calculate the leading nuclear structure factors that exhibit the coherent contribution of several nucleons in the nucleus. We call these generalized SI interactions. We consider both WIMP couplings to one and two nucleons, with special emphasis on the two-nucleon contributions related to the diagonal and non-diagonal parts of the energy-momentum tensor. We note that at this point our structure factors are not fully consistent, in the sense that our many-body calculations are based on phenomenological interactions instead of ChEFT. Such consistent studies are presently only available for few-nucleon systems Korber:2017ery . In addition, nuclear response functions based on nuclear states obtained using ChEFT interactions are available for light nuclei Gazda:2016mrp . However, such very light isotopes are not used in leading direct detection experiments. While Ref. Hoferichter:2016nvd was limited to a xenon target, here we study all the stable isotopes of fluorine, silicon, argon, and germanium, the nuclear targets used and considered in present and future direct detection experiments. In addition, we provide a Python notebook to facilitate the use of our structure factors for both theorists and experimentalists.

The rest of the article is organized as follows. In Sec. II we describe the formalism and discuss which terms we include in the cross section, focusing on a spin-$1/2$ WIMP. Section III introduces the nuclear structure many-body calculations that describe the nuclear targets used in direct-detection experiments. The results for our calculated structure factors are given in Sec. IV. We distinguish between WIMP couplings to one nucleon, discussed in Sec. IV.1, and the coupling to two nucleons, the subject of Sec. IV.2. The most important new aspects of our implementation of two-body effects are highlighted in Sec. V. We conclude with a summary of the main findings of this work in Sec. VI. Details of the matching to NREFT, the nucleon matrix elements, matching for a spin-[math] WIMP, and the calculation of two-body interactions are provided in Appendices A–C. All our results are also available in the form of a Python notebook, with a brief User’s guide in Appendix D.

II Formalism

Based on Ref. Hoferichter:2016nvd we consider the following cross section for the generalized SI WIMP–nucleus scattering:

[TABLE]

where $q=|\mathbf{q}|$ is the three-momentum transfer between the WIMP and the nucleus, and $m_{\chi}$, $m_{N}$, and $m_{\mathcal{N}}$ denote the WIMP, nucleon, and nuclear masses, respectively. The relative velocity of the WIMP is $v=|\mathbf{v}|$, while $(v^{\perp}_{T})^{2}=v^{2}-q^{2}/(4\mu_{\mathcal{N}}^{2})$ Fitzpatrick:2012ix is the velocity of the WIMP with respect to the nucleus, with reduced mass $\mu_{\mathcal{N}}=m_{\mathcal{N}}m_{\chi}/(m_{\mathcal{N}}+m_{\chi})$. Each term in Eq. (II) contains a structure factor $\mathcal{F}$ and a coupling $c$. The structure factors $\mathcal{F}$ encode the nuclear structure aspects of the scattering, and are the main subject of this work. The couplings $c$ are a convolution of the Wilson coefficients, which describe the fundamental interaction of the WIMPs with the quarks and gluons, and the hadronic matrix elements. They therefore depend on the particular scenario for physics beyond the Standard Model considered. Finally, the kinematic factors

[TABLE]

with the reduced mass $\mu_{N}=m_{N}m_{\chi}/(m_{N}+m_{\chi})$, appear in the contribution from the subleading NREFT operators $\mathcal{O}_{5,8,11}$, defined in Appendix A.

II.1 Effective Lagrangian and couplings

To be definite, we consider the case of a spin-$1/2$, Standard-Model singlet $\chi$ throughout the main text, but to demonstrate that the decomposition in Eq. (II) applies in full generality we also provide the matching relations for a spin-[math] WIMP, see Appendix B. For spin-$1/2$, the relevant terms in the effective Lagrangian Goodman:2010ku ; Drees:1993bu ; Belanger:2008sj are

[TABLE]

where we have only listed operators that can lead to coherently enhanced responses or feature in the standard SD interaction, and $C_{q}^{VV}=C_{q}^{TT}=\tilde{C}_{q}^{TT}=0$ for a Majorana particle. In $\mathcal{L}_{\chi}^{(7)}$ we already integrated out the heavy quarks Shifman:1978zn , whose effect is absorbed in

[TABLE]

while elsewhere the sum runs, in principle, over all quark flavors $q$. $C^{S}_{g}$ and $\tilde{C}^{S}_{g}$ are the original coefficients of the gluon operator $\bar{\chi}\chi\,\alpha_{s}G_{\mu\nu}^{a}G^{\mu\nu}_{a}$ and $\bar{\chi}i\gamma_{5}\chi\,\alpha_{s}G_{\mu\nu}^{a}G^{\mu\nu}_{a}$, respectively, rewritten in terms of the trace of the energy-momentum tensor:

[TABLE]

Moreover, we introduced its traceless components $\bar{\theta}^{\mu\nu}=\bar{\theta}^{\mu\nu}_{q}+\bar{\theta}^{\mu\nu}_{g}$ in the context of the dimension-$8$ spin-$2$ contribution:

[TABLE]

with covariant derivative $D_{-}^{\mu}=\overrightarrow{D}^{\mu}-\overleftarrow{D}^{\mu}$ and symmetrizer $A^{\{\mu}B^{\nu\}}=(A^{\mu}B^{\nu}+A^{\nu}B^{\mu})/2$.

In detail, the coefficients in Eq. (II) are given by (nucleons $N$ encompass protons and neutrons, $N=p$ or $n$)

[TABLE]

where $\zeta=1(2)$ for a Dirac (Majorana) particle.111The matching relations thus include the Majorana symmetry factor $2$. In the Dirac case, the relative sign between scalar and vector coefficients changes for the cross section of $\bar{\chi}$–nucleon scattering Belanger:2008sj .

The various couplings refer to combinations of Wilson coefficients from Eq. (3) and nucleon matrix elements in the scalar ($f_{N}$, $\dot{f}_{N}$, $\tilde{f}_{N}$), vector ($f_{1}^{V,N}$, $\dot{f}_{1}^{V,N}$, $f_{2}^{V,N}$, $\tilde{f}_{1}^{V,N}$), tensor ($f^{T,N}$, $\tilde{f}^{T,N}$), and spin-$2$ ($f_{N}^{(2)}$) channels. The explicit expressions, together with the matching onto NREFT, are collected in Appendix A.

II.2 Pion matrix elements

Here, we highlight the couplings to the pion in Eq. (II.1):

[TABLE]

In this case, the relevant matrix elements are the following: first, the scalar couplings to $u$- and $d$-quarks are given by

[TABLE]

where $\xi_{ud}=(m_{d}-m_{u})/(m_{d}+m_{u})=0.37(3)$ with $m_{u}/m_{d}=0.46(3)$ from Ref. Aoki:2016frl (and $f_{s}^{\pi}\approx 0$). Second, the coupling to the trace anomaly is only sensitive to the sum $f_{u}^{\pi}+f_{d}^{\pi}=1$, so that no new coupling appears in $f_{\pi}^{\theta}$. The spin-$2$ couplings are related to moments of pion parton distribution functions (PDFs) $q_{\pi}(x)$:

[TABLE]

and similarly for the gluon, subject to the constraint

[TABLE]

The pion PDFs have often been studied assuming a single valence and sea distribution with moments $\langle x\rangle_{\text{v}}^{\pi}$ and $\langle x\rangle_{\text{s}}^{\pi}$ Badier:1983mj ; Sutton:1991ay ; Wijesooriya:2005ir , in such a way that

[TABLE]

and by means of the sum rule in Eq. (11)

[TABLE]

With $f^{(2)}_{g,\pi}=0.47(15)$ Badier:1983mj and $\langle x\rangle_{\text{v}}^{\pi}=0.217(11)$ Wijesooriya:2005ir this gives

[TABLE]

consistent with a recent calculation in lattice QCD $\langle x\rangle_{\text{v}}^{\pi}=0.214(15)(^{+12}_{-9})$ Abdel-Rehim:2015owa .222Note that all spin-$2$ couplings are scale-dependent quantities, so that, at a scale of $\mu=2\,\text{GeV}$ in $\overline{\text{MS}}$, the value from Ref. Wijesooriya:2005ir actually becomes $\langle x\rangle_{\text{v}}^{\pi}=0.256(13)$ Abdel-Rehim:2015owa . This has been taken into account in Eq. (14). A more recent, global analysis finds at $\mu=2\,\text{GeV}$ Barry:2018ort

[TABLE]

in agreement with Eq. (14), but considerably more precise. For completeness, we reproduce the analogous expressions for the nucleon matrix elements in Appendix A.

II.3 Dipole operators

For a Dirac WIMP, a possible extension beyond Eq. (3) concerns dark matter candidates with a non-vanishing dipole moment, corresponding to the effective dimension-$5$ Lagrangian Barger:2010gv ; Banks:2010eh

[TABLE]

where due to $\sigma^{\mu\nu}\gamma_{5}=\frac{i}{2}\epsilon^{\mu\nu\alpha\beta}\sigma_{\alpha\beta}$ the second term involving the dual field strength tensor $\tilde{F}^{\mu\nu}=\frac{1}{2}\epsilon^{\mu\nu\alpha\beta}F_{\alpha\beta}$ is equivalent to $-\bar{\chi}\sigma^{\mu\nu}i\gamma_{5}\chi F_{\mu\nu}$. These operators produce long-range tree-level interactions via the exchange of a photon, which leads to the NR one-body amplitudes

[TABLE]

Here, $F_{1/2}^{N}(t)$ are the Dirac/Pauli form factors of the nucleon with the full dependence on the relativistic momentum transfer $t=-q^{2}$ (up to relativistic corrections). The corresponding extension of Eq. (II) is straightforward, with the photon poles introducing new terms, besides additional contributions to some of the existing ones. Formally, this amounts to ($q$-dependent) terms in the matching relations in Eq. (II.1):

[TABLE]

where we have already taken $\zeta=1$ due to the absence of tensor currents for Majorana particles. The Dirac-form-factor radii are related to the Sachs ones by

[TABLE]

In the spirit of the present study, Eq. (18) neglects the non-coherent contributions from $\mathcal{O}_{4}$ and $\mathcal{O}_{6}$. Further, there are many more dimension-$7$ operators involving the (electroweak) field strength tensors Crivellin:2014gpa ; Brod:2017bsw , and matching relations similar to Eqs. (II.3) and (18) could be extended accordingly.

In the form (II.3) the amplitudes automatically include radius corrections, subsumed in the full electromagnetic form factors. In addition, Eq. (16) could in principle produce new two-body currents. In the end, such terms take a similar form as the axial-vector–vector two-body currents identified in Ref. Hoferichter:2015ipa , and are only suppressed by two chiral orders compared to the leading pieces of Eq. (II.3). However, as argued in Ref. Hoferichter:2016nvd , after summation over spins their isospin structure $[\boldsymbol{\tau}_{1}\times\boldsymbol{\tau}_{2}]^{3}$ leaves only an isovector coherent enhancement suppressed by $(N-Z)/A$ with respect to the scalar two-body current. We will continue to neglect such effects in the remainder of this work.

II.4 Scaling of operators

The global picture that arises in this way is summarized in Table 1. For each nuclear structure factor we list: the relativistic operators that contribute, the corresponding NREFT operator if applicable, the chiral scaling as well as the overall scaling including coherence, and finally whether or not the respective contribution interferes with the leading $\mathcal{O}_{1}$ operator.

Table 1 reflects the strategy underlying the present work, in that it includes all contributions, either two-body currents or subleading one-body operators, that appear up to third order in the chiral expansion and receive some form of coherent enhancement. Selected coherent chiral fourth-order contributions are shown as well, mainly because in many cases this is the order when some of the coherent NREFT operators first enter. For comparison, the table also includes the leading relativistic operators that produce the NR expansion related to the standard SD interactions, even if these are not coherent.

III Nuclear structure calculations

The calculation of the nuclear structure factors requires a many-body approach that describes the ground states of the target nuclei considered. As in previous works Menendez:2012tm ; Klos:2013rwa ; Baudis:2013bba ; Vietze:2014vsa ; Hoferichter:2016nvd we use the nuclear shell model, one of the most successful many-body approaches in medium-mass and heavy nuclei Caurier:2004gf . For all calculations we have used the shell model code ANTOINE Caurier:1999 ; Caurier:2004gf .

The shell model is based on the solution of the quantum many-body problem in a reduced configuration space where the Schrödinger equation for the nuclear ground and low-energy excited states can be solved exactly. We highlight two important aspects. First, the configuration space used in the calculation needs to capture the nuclear structure properties relevant for the process of interest. The limitation of using a restricted configuration space stems from the difficulty to solve the nuclear many-body problem in a nontruncated space for heavier nuclei. For nuclear targets used in direct detection experiments, calculations of structure factors without such truncations exist up to 4He Korber:2017ery ; Gazda:2016mrp ; Andreoli:2018etf and could be performed in the near future up to 40Ca. Second, calculations must use an effective interaction appropriate for such a configuration space.

Until very recently, nuclear shell-model calculations relied on phenomenological effective interactions. In spite of being derived from nucleon–nucleon ($NN$) scattering data, these interactions have to be adjusted phenomenologically, mostly the part that describes the single-particle aspects of the nuclear interaction, referred to as monopole part, in order to achieve a better agreement with the nuclear structure of heavier nuclei. Progress in nuclear theory has improved this picture by including in the starting point of the derivation of effective interactions, in addition to $NN$ interactions, also nuclear forces between three nucleons, 3$N$ forces. Three-nucleon forces are the analog of two-body currents in the coupling to external probes, such as WIMPs, to nucleons. Consistent $NN$ and 3$N$ interactions can be derived in the ChEFT framework Hammer:2012id ; Machleidt:2011zz ; Epelbaum:2008ga . Starting from these ChEFT $NN$ and 3$N$ forces it is possible to derive effective interactions that, without further phenomenological adjustments, reproduce nuclear spectroscopy rather well in nuclei with nucleon number $A$ comparable to direct-detection nuclear targets Hebeler:2015hla ; Hergert:2015awm ; Hagen:2013nca . First studies have just started to extend these techniques to study electromagnetic and weak transitions in medium-mass nuclei Hagen:2015yea ; Parzuchowski:2017wcq ; Morris:2017vxi . In addition, the effective-theory character of ChEFT, combined with the consistency between nuclear forces and currents, provides a framework to quantify nuclear-structure uncertainties Epelbaum:2014efa ; Furnstahl:2015rha ; Carlsson:2015vda . In this work, however, we follow the standard approach and use phenomenological effective interactions. This prevents the assessment of reliable nuclear-structure uncertainties, which will become possible in the future.

The configuration spaces and effective interactions we have used are described below. All single-particle orbitals belong to a three-dimensional harmonic oscillator basis $nl_{j}$, where $n$ is the principal quantum number, and $l,j$ denote the orbital and total angular momentum.

The lightest target we have studied has only one stable isotope, 19F, which we already considered in Ref. Klos:2013rwa . Here we use the same USDB effective interaction and the $0d_{5/2}$, $1s_{1/2}$, and $0d_{3/2}$ single-particle orbitals, for both neutrons and protons. This configuration space is known as the $sd$ shell. Reference Klos:2013rwa used the same configuration space and effective interaction for 29Si, and here we extend the study to the other two stable silicon isotopes, 28,30Si.

For argon, we study only 40Ar. The configuration space we consider is significantly larger, and comprises seven single-particle orbitals for neutrons and protons, the three $sd$-shell orbitals plus the $0f_{7/2}$, $1p_{3/2}$, $0f_{5/2}$, and $1p_{3/2}$, where the latter comprise the so-called $pf$ shell. We use the SDPF.SM effective interaction, which describes well the electromagnetic properties of ground states and the coexistence of spherical and deformed states in this mass region Caurier:2007ee ; Ruiz:2015yra . In order to make the diagonalizations in the configuration space feasible, we need to truncate our many-body calculations, by keeping the $0d_{5/2}$ orbital filled with nucleons, and restricting the number of excitations from $sd$-shell to $pf$-shell orbitals to 8. These are similar truncations to those in Ref. Caurier:2007ee , and limit the dimension of the diagonalization to $5\cdot 10^{8}$. We have also preformed calculations of the stable isotopes 36,38Ar which are of the same quality as those for 40Ar. Nonetheless we have not included them in our study, because the natural abundance of these isotopes is very minor with less than $0.3\%$.

Finally, there are five stable stable germanium isotopes 70,72,73,74,76Ge. Consistently with our study of 73Ge in Ref. Klos:2013rwa we use the RG effective interaction Menendez:2009xa in a configuration space consisting on the $1p_{3/2}$, $0f_{5/2}$, $1p_{3/2}$, and $g_{9/2}$ single-particle orbitals. We have performed calculations with alternative effective interactions in the same configuration space Caurier:2007wq ; Menendez:2008jp ; Honma:2009zz ; jj44b , and while the excitation spectra may be somewhat different to the ones predicted by the RG interaction, the impact on the nuclear structure factors is very small at the momentum transfers relevant to direct detection searches.

Figures 1–6 compare the low-energy excitation spectra of the stable isotopes of fluorine, silicon, argon, and germanium with our theoretical predictions. In all cases, our calculations are in very good agreement with experiment, especially for nuclei with even number of nucleons. In some cases, especially for the odd-mass nucleus 73Ge, some experimental states are not well reproduced. This is very likely due to the limitation of the configuration space used in our calculation, because the description of 73Ge is of similar quality with the other effective interactions we have studied.

To complement the assessment of the quality of the nuclear structure calculations, Tables 2 and 3 compare theoretical and experimental electromagnetic observables for all these nuclear targets. The comparison includes nuclear charge radii, electromagnetic moments for ground and lowest-excited states, and nuclear matrix elements for selected electromagnetic transitions between low-lying states. Details on the shell model calculation of nuclear moments and matrix elements can be found, e.g., in Ref. Caurier:2004gf . Charge radii Bertozzi:1972jff , the properties most relevant for coherent nuclear structure factors, are very well reproduced, better than $3\%$ in all cases, and better than $1\%$ in the heavier argon and germanium. For fluorine and silicon, Table 2 shows an excellent agreement between theory and experiment, as the majority of the predictions reproduce data within experimental uncertainties. For argon and germanium, Table 3 also shows a reasonable agreement of the nuclear structure calculations with experiment. Nuclear matrix elements within the same isotope can vary over two orders of magnitude, and the theoretical results reproduce well the corresponding hierarchy. In germanium, some theoretical electric moments and transitions underestimate experiment moderately. This suggests that the configuration space used in the calculation is not sufficient to fully account for the most collective states, consistently with the findings from the comparison of the 73Ge excitation spectrum. We do not expect, however, any significant effect in the ground states involved in elastic WIMP–nucleus scattering.

In Sec. IV and the Python notebook we also cover nuclear structure factors for xenon, which we studied in the context of coherent SI scattering in Refs. Vietze:2014vsa ; Hoferichter:2016nvd , and more generally in Refs. Klos:2013rwa ; Menendez:2012tm . The quality of the nuclear structure calculations of stable xenon isotopes is similar to that of argon or germanium, the heaviest nuclear targets considered in this section. Calculations for xenon have been compared to experimental data in Refs. Vietze:2014vsa ; Klos:2013rwa for excitation spectra, and in Refs. Menendez:2012tm ; Sieja:2009zz for electromagnetic properties. Given that the charge radii are closely connected to the nuclear structure factors considered here, Table 4 compares theoretical results with experiment. In all isotopes the calculations reproduce measured radii to better than 1%.

IV Structure factors

Apart from the dependence on $q$, $m_{N}$, and $m_{\chi}$ that is predicted by ChEFT, the generalized SI WIMP–nucleus cross section in Eq. (II) depends on six independent structure factors $\mathcal{F}$. Four of them correspond to the coupling of the WIMP to one nucleon, which can be the same for protons and neutrons (as in the two isoscalar structure factors) or opposite (as in the two isovector ones). In addition, two independent structure factors characterize the simultaneous coupling of WIMPs to two nucleons. In this section we evaluate these six structure factors for the nuclear targets considered in present and future direct detection experiments. An overview of the various contributions, excluding interference terms, is provided in Fig. 7. In the following, we present in detail our results for the one-body (1b) and two-body (2b) structure factors.

IV.1 One-body structure factors

As discussed in Refs. Fitzpatrick:2012ix ; Anand:2013yka there are two different nuclear responses describing the coupling of a WIMP to a single nucleon, $\mathcal{F}^{M}$ and $\mathcal{F}^{\Phi^{\prime\prime}}$, which receive coherent contributions from several nucleons in the nucleus. In addition to the standard SI scattering, the nuclear response $\mathcal{F}^{M}$ describes a subleading contribution that corresponds to the NREFT operators $\mathcal{O}_{5,8,11}$, but by far the most important among them is $\mathcal{O}_{11}$. In addition, the so-called radius correction to the standard SI structure factor is also coherent Hoferichter:2016nvd . Dropping the contributions from $\mathcal{O}_{5,8}$, the scattering cross section including one-nucleon couplings simplifies to

[TABLE]

where $\tilde{c}_{I}^{M}$ can be identified with $c_{I}^{M,11}$ from Eq. (II). Even though there are only four independent isoscalar contributions (plus four isovector ones), in the most general case where all contributions in Eq. (IV.1) compete, the interference of all of them generates a plethora of individual terms that could be considered.

Figures 8, 9, 10, and 11 show the leading contributions to the cross section for 132Xe, 74Ge, 40Ar, and 19F, respectively. The results for xenon use the results of Ref. Hoferichter:2016nvd . Figures 8, 9, 10, and  11 assume that all couplings $c$ are equal to $1$, and for the NREFT $\mathcal{O}_{11}$ term $m_{\chi}=2$ GeV, which implies that for heavier WIMPs the importance of this term will always be smaller than in the figures.

Figures 8, 9, 10, and 11 highlight that the standard SI contribution (solid black lines) is indeed expected to be dominant. Moreover, leading corrections to generalized SI scattering come from the interference of the standard SI and other terms, such as its isovector counterpart (dotted-dashed black), the isoscalar radius corrections (solid blue), or the isoscalar $\mathcal{F}^{\Phi^{\prime\prime}}$ term (solid red). The only exceptions are, first, the purely isovector SI structure factor (dashed black), suppressed by one or two orders of magnitude. Second, the $\mathcal{O}_{11}$ contributions (solid green), suppressed by around four orders of magnitude in all cases (the suppression will be larger for heavy WIMPs $m_{\chi}>2\,\text{GeV}$).

The variable importance of these contributions is set by the nuclear structure of the corresponding nuclear targets. Isovector contributions are relatively more important in neutron-rich xenon than in the $N\approx Z$ fluorine. On the other hand, the $\mathcal{F}^{\Phi^{\prime\prime}}$ contributions are relatively more important in the heavier xenon and germanium, because these targets have more nucleons in single-particle orbitals with aligned spin and orbital angular momentum. In contrast, the $\mathcal{F}^{\Phi^{\prime\prime}}$ contributions are more suppressed in lighter targets such as argon and especially fluorine, which tend to have nucleons more equally distributed in orbitals with spin parallel and antiparallel to the orbital angular momentum.

IV.2 Two-body structure factors

The two-body amplitudes of the three scalar channels, scalar–scalar ($SS$), trace anomaly $\theta^{\mu}_{\mu}$ (“$\theta$” in short), and spin-$2$ (referred to by “$(2)$”), read

[TABLE]

where the couplings $f_{\pi}$, $f_{\pi}^{\theta}$, $f_{\pi}^{(2)}$ are defined in Eq. (II.2). $\boldsymbol{\sigma}_{i}$ and $\boldsymbol{\tau}_{i}$ refer to spin and isospin operators for nucleon $i$, $F_{\pi}$ is the pion decay constant, and $g_{A}$ the axial charge of the nucleon. Throughout, we use PDG values Tanabashi:2018oca , except for the particle masses, for which we use isospin averages of $m_{N}=938.92\,\text{MeV}$ and $M_{\pi}=138\,\text{MeV}$.

While Refs. Cirigliano:2012pq ; Hoferichter:2015ipa ; Hoferichter:2016nvd introduced coherent two-body currents, the present work includes for the first time the contact-term contributions to the $\theta^{\mu}_{\mu}$ response involving $C_{S}$ and $C_{T}$ and the entire spin-$2$ contribution. In particular the consistent inclusion of the contact operators is a crucial improvement in this work (the result for $\theta^{\mu}_{\mu}$ was already used in Ref. Hoferichter:2017olk ), see Sec. V for an extended discussion.

By including the relativistic corrections of subleading one-body terms, both in the $\theta^{\mu}_{\mu}$ and spin-$2$ channels, it is possible to write the three physical responses in terms of just two new structure factors: see Sec. V and Appendix C for more details and the precise definition of contact operators. The two structure factors in the naive (non-interacting) shell model read

[TABLE]

where the second structure factor is normalized as $\mathcal{F}_{\text{b}}(0)=-2E_{\text{b}}/M_{\pi}$ with the binding energy of the nucleus $E_{\text{b}}<0$. $\mathcal{F}_{\pi}(q^{2})$ corresponds to the scalar–scalar two-body current, and by defining $\mathcal{F}_{\text{b}}(q^{2})$ as above, the three physical channels ($SS$, $\theta$, $(2)$) are described in terms of these two structure factors because

[TABLE]

Figures 12, 13, 14, and 15 show the structure factors for 132Xe, 74Ge, 40Ar, and 19F that include contributions from the coupling to two nucleons. Scalar couplings are described by the $\mathcal{F}_{\pi}$ contributions (thick solid orange line), which can interfere with the SI contribution (thin solid orange) and its isovector counterpart (dotted-dashed orange) and with an independent two-nucleon coupling (maroon dotted-double-dashed).

The two-nucleon coupling to the trace anomaly receives two contributions. According to Eq. (IV.2), the first one can be described by $\mathcal{F}_{\pi}$ and the second one by the structure factor $\mathcal{F}_{b}$. Figures 12, 13, 14, and 15 show the $\mathcal{F}_{\pi}$, $\mathcal{F}_{b}$ structure factors, and their interferences with one-nucleon couplings. In particular, besides the terms described above, the figures show the full $\mathcal{F}_{b}$ structure factor (thick solid indigo line), its interference with the SI term (thin solid indigo), and its isovector counterpart (dotted-dashed indigo), and finally the interference with the radius correction term (dotted indigo).

Similarly, spin-$2$ two-nucleon couplings contribute via $\mathcal{F}_{\pi}$ and $\mathcal{F}_{b}$ terms according to Eq. (IV.2). Figures 16 and 17 compare the contribution of the physical combination of the $\mathcal{F}_{\pi}$ and $\mathcal{F}_{b}$ structure factors that originate in the scalar, trace anomaly, and spin-$2$ two-nucleon couplings. Taking all coefficients $c$ to unity, the dominant effect is given by the trace anomaly (dashed line), followed by the spin-$2$ term (dotted-dashed), and the scalar contribution (solid). However, we stress that this hierarchy only reflects the nuclear structure aspects and does not need to be followed by particular models with definite Wilson coefficients. The nucleon matrix elements and especially the BSM couplings can alter the importance of the different structure factors in Figs. 16 and 17.

V Contact operators in two-body structure factors

The results for the two-body currents presented in Sec. IV.2 are based on the ChEFT formalism developed in Ref. Hoferichter:2016nvd , see Appendix C for more details. In particular, the chiral power counting follows the proposal from Ref. Weinberg:1990rz ; Weinberg:1991um , in which the scaling of operators is estimated by dimensional analysis. In this way, the leading contribution for a scalar current stems from pion-exchange diagrams, since $(N^{\dagger}N)^{2}$ contact operators require an additional insertion of a scalar source that is counted in the same way as the quark mass matrix and therefore only appears at subleading order, see, e.g., Ref. Epelbaum:2002gb . In alternative formulations Kaplan:1998tg ; Kaplan:1998we where contact operators are promoted to lower orders, the pion-exchange diagrams would be accompanied by additional contact-term contributions at the same order, as also suggested by renormalization group arguments for external currents Valderrama:2014vra . A conclusive test of the importance of these contact operators would require detailed studies of the scalar current in light nuclei along the lines of Ref. Korber:2017ery , using recent precision chiral potentials Epelbaum:2014sza ; Entem:2017gor ; Reinert:2017usi , contrasted to nuclear $\sigma$-terms from lattice QCD Beane:2013kca ; Chang:2017eiq . Work along these lines is in progress.

In practice, the question arises as to how to deal with such potential contact operators in nuclear many-body calculations. In fact, even in the Weinberg power counting contact operators do occur at leading order in the coupling to the trace anomaly $\theta^{\mu}_{\mu}$ and a spin-$2$ source. These terms were neglected in Ref. Hoferichter:2016nvd , so that the resulting $\theta^{\mu}_{\mu}$ response was incomplete (their contributions were, however, included in Ref. Hoferichter:2017olk ). Here, we show how in the Weinberg power counting the contact operators in these channels are canonically renormalized in terms of nuclear binding energies, a mechanism that no longer applies once additional contact terms are introduced.

V.1 Trace anomaly

The simplest example for the chiral realization of the energy-momentum tensor $\theta_{\mu\nu}$ in ChEFT is the tree-level pion matrix element:

[TABLE]

i.e., for an on-shell pion

[TABLE]

One-loop corrections have been worked out in Ref. Donoghue:1991qv . Likewise the diagrams in Fig. 18 for the nucleon give

[TABLE]

with loop integral

[TABLE]

In the two-nucleon sector, see Fig. 19, the first term in Eq. (22) follows from the pion-exchange diagram $(a)$ by means of Eq. (26) (diagram $(b)$ only enters at higher orders). An additional contribution arises from the contact-term Lagrangian:

[TABLE]

where the second term derives from the NR expansion of the spin vector $S^{\mu}=(0,\boldsymbol{\sigma}/2)$. The corresponding term in the trace

[TABLE]

gives the second term in Eq. (22), represented by diagram $(c)$ in Fig. 19.

In contrast to the scalar current, the pion-exchange piece in Eq. (22) behaves as a contact term for $q_{i}\to\infty$, due to the momentum dependence of the pion coupling, see Eq. (26). As first noted in Ref. Hoferichter:2017olk , for vanishing momentum transfer these two pieces combine to the $NN$ potential $V_{NN}$. Together with the kinetic-energy operator $T$, this suggests the renormalization prescription

[TABLE]

where $E_{\text{b}}<0$ is the binding energy of the nucleus, whose wave function $|\Psi\rangle$ should be obtained from $NN$ interactions only. In practice, the kinetic-energy operator formally enters at higher orders but arises naturally from relativistic corrections, while $3N$ forces correspond to higher-order corrections. In order to determine $C_{S}$ and $C_{T}$ we therefore use experimental binding energies, corrected for Coulomb interactions according to Ref. Duflo:1995ep . Details of the renormalization both for $\theta^{\mu}_{\mu}$ and spin-$2$ terms are given in Appendix C.

For the numerical analysis we ignore $C_{T}$ contributions, because the magnitude of $C_{T}$ is expected to be small due to approximate $SU(4)$ symmetry of $NN$ interactions at low energies, with corrections that can be shown to be suppressed by $1/N_{c}^{2}$ Kaplan:1995yg ; Kaplan:1996rk ; Mehen:1999qs . The resulting $C_{S}$ values in Table 5 are consistent with the expectation from dimensional analysis Epelbaum:2004fk :

[TABLE]

where we have used that $C_{T}=0$ implies $\tilde{C}_{{}^{1}S_{0}}=\tilde{C}_{{}^{3}S_{1}}$. Moreover, the values in Table 5 also agree with typical fits to the $NN$ system, e.g., at LO Ref. Epelbaum:2014efa finds $C_{S}=(-56.5\ldots-118.3)\,\text{GeV}^{-2}$ for cutoffs in the range $R=(0.8\ldots 1.2)\,\text{fm}$.

V.2 Spin 2

Despite only entering at dimension-$8$ in Eq. (3), spin-$2$ contributions become relevant, for instance, in the context of heavy WIMPs, where significant cancellations with spin-[math] terms have been observed Hill:2013hoa , enhancing the importance of higher-order corrections. The relevant operators are the traceless parts of the energy-momentum tensor, given in Eq. (II.1). As can be seen from the WIMP part of ${\mathcal{L}}^{(8)}_{\chi}$ in Eq. (3), the dominant contribution arises from the $\mu=\nu=0$ components, leading to the amplitudes $\mathcal{M}_{1,\text{NR}}^{(2)}$ in Eq. (A) and $\mathcal{M}_{2,\text{NR}}^{(2)}$ in Eq. (23). An extension to subleading components is straightforward, since due to Lorentz invariance the pion-exchange contribution becomes proportional to $q_{\mu}q_{\nu}-\frac{g_{\mu\nu}}{4}q^{2}$. Then the full expression can be reconstructed from the $00$ component, which we have identified with $\mathcal{F}^{(2)}_{\pi}(q^{2})$.

The chiral realizations of the matrix elements of $\bar{\theta}^{\mu\nu}_{q,g}$ have been studied in detail in the literature, both for the pion and the nucleon Arndt:2001ye ; Chen:2001nb ; Detmold:2005pt ; Diehl:2005rn ; Ando:2006sk ; Diehl:2006js ; Wein:2014wma . Here, we only retain the leading couplings related to moments of pion and nucleon PDFs, resulting in the one- and two-body contributions $\mathcal{M}_{1,\text{NR}}^{(2)}$ and $\mathcal{M}_{2,\text{NR}}^{(2)}$, as well as the relativistic corrections in Eq. (55). Motivated by the EMC effect Aubert:1983xm , similar methods have been applied in the context of spin-$2$ couplings in multi-nucleon systems Chen:2004zx ; Chen:2016bde . Therefore, measurements of nuclear PDFs could, in principle, provide independent cross checks on the resulting spin-$2$ structure factor. However, in practice this is not possible with currently employed parameterizations Hirai:2007sx ; deFlorian:2011fp ; Kovarik:2015cma ; Eskola:2016oht : in the dark matter context, the main effect of the two-body corrections modifies the normalization at $q=0$ away from the fully coherent single-particle expectation $\mathcal{F}(0)=A$, see Cirigliano:2012pq ; Hoferichter:2016nvd for the scalar channel. Presently, nuclear PDFs $q^{A}(x)$ are studied based on bound-proton PDFs $q^{p/A}(x)$ restricted onto the range $x\in[0,1]$, in such a way that the full PDF is reconstructed by $q^{A}(x)=Zq^{p/A}(x)+Nq^{n/A}(x)$. The moments of this $q^{A}(x)$ are therefore, by definition, normalized to the coherent limit $A$ and cannot be used to cross check the two-body effect that we derived from the spin-$2$ couplings of the pion.

VI Summary

We have presented a comprehensive analysis of the generalized SI scattering of spin-$1/2$ and spin-[math] WIMPs off atomic nuclei. Our analysis considers all contributions that can receive the coherent enhancement from several nucleons in the nucleus, keeping terms up to third order in ChEFT. This includes both the coupling of WIMPs to one nucleon as well as to two nucleons. For two-body interactions we provide, for the first time, a full and consistent treatment of the contact operators that appear at the same order as pion-exchange diagrams, arguing that these contributions can be renormalized to the nuclear binding energy. As a result, just two nuclear structure factors are enough to characterize two-body interactions via scalar operators, the trace of the energy-momentum tensor, and spin-$2$ operators.

Taking into account all these contributions, we give all one-body and two-body nuclear structure factors relevant for the coherent WIMP scattering off fluorine, silicon, argon, germanium and xenon, covering the targets of the most advanced direct detection searches. For that purpose, we perform large-scale nuclear shell model calculations with configuration spaces and nuclear interactions that describe very well the structure of these nuclei.

Our analysis identifies the parameters that can, at least in principle, be separately constrained in direct detection experiments. These parameters subsume both the BSM couplings of WIMPs with quarks and gluons and the hadronic matrix elements that embed these quark-level operators into hadrons. The corresponding matching relations are illustrated in detail for both spin-1/2 and spin-0 WIMPs.

The main results of our work, encoded in the nuclear structure factors and the relation between direct-detection experiments and BSM couplings, are available as supplementary material in a Python notebook. These results form the basis for a comprehensive study of WIMP–nucleus interactions based on ChEFT. Future extensions concern non-coherent WIMP–nucleus interactions, for which more parameters and nuclear structure factors need to be considered. Accordingly, if the coherent contributions studied in this paper are strongly suppressed, the identification of the underlying quark-level interactions becomes even more challenging. On the other hand, progress in ab initio nuclear theory paves the way towards fully consistent structure factors from many-body calculations based on ChEFT Andreoli:2018etf ; Hebeler:2015hla ; Hergert:2015awm ; Hagen:2013nca ; Korber:2017ery ; Gazda:2016mrp . Such improved nuclear structure factors, including their momentum-dependence, will further help distinguish among possible BSM scenarios.

Acknowledgements.

We thank V. C. Antochi, P. Barry, W. Detmold, E. Epelbaum, A. Fieguth, H.-W. Hammer, C. Hasterok, D. B. Kaplan, K. Kovařík, H. Krebs, T. Marrodan, F. Olness, and J. de Vries for valuable discussions. This work was supported in part by the US DOE (Grant No. DE-FG02-00ER41132), the National Science Foundation (Grant No. NSF PHY-1748958), the ERC (Grant No. 307986 STRONGINT), the DFG through SFB 1245 (Projektnummer 279384907), the Max-Planck Society, the Japanese Society for the Promotion of Science KAKENHI through grant 18K03639, MEXT as “Priority Issue on Post-K computer” (Elucidation of the fundamental laws and evolution of the universe), JICFuS, and the CNS-RIKEN joint project for large-scale nuclear structure calculations. J. M. and A. S. thank the Institute for Nuclear Theory at the University of Washington for its hospitality and the US DOE for partial support.

Appendix A Matching to NREFT and nucleon matrix elements

For the matching onto NR single-nucleon operators we use the conventions

[TABLE]

with

[TABLE]

and

[TABLE]

as well as the operator basis

[TABLE]

where $q=|\mathbf{q}|$.

The expressions above refer to the WIMP–nucleon system. In the nucleus, the operator $\mathbf{v}^{\perp}$ generates two kinds of contributions Fitzpatrick:2012ix ; Anand:2013yka . First, there are operators dependent on the WIMP velocity with respect to the nucleus, $\mathbf{v}^{\perp}_{T}$ [see the terms involving $\xi_{i}(q,v^{\perp}_{T})$ in Eq. (II)]. These terms are very suppressed in the scattering amplitude because $v^{\perp}_{T}=|\mathbf{v}^{\perp}_{T}|\sim 10^{-3}$. On the other hand, $\mathbf{v}^{\perp}$ generates terms that contain the nucleon velocity operator. In this case, the operators are mildly suppressed by $q/m_{N}$, and include an additional derivative. The latter terms are fully responsible for the $\mathcal{F}^{\Phi^{\prime\prime}}$ structure factor.

Back to the WIMP–nucleon level, the coherently enhanced terms listed in Table 1, complemented by the leading SD response, are derived from

[TABLE]

where we have ignored the non-coherent terms in the NREFT expansion. For the remaining couplings the momentum dependence is indicated by the relativistic momentum transfer $t$, which reduces to $t=-q^{2}$ up to relativistic corrections.

Some of the amplitudes in Eq. (A) receive contributions that break Galilean invariance. For such terms, which only appear beyond $\mathcal{O}(p^{3})$ in the chiral expansion, we assume center-of-mass kinematics, for which the velocity in Eq. (36) simplifies to

[TABLE]

In principle, corrections to this identification would need to be considered when calculating the nuclear structure factors, similar to the boost correction in Ref. Beane:2002wk , but given that these contributions are already highly suppressed, we only keep the center-of-mass component. The appearance of such Galilean-invariance-breaking terms at subleading orders in the NR expansion has been pointed out in Ref. Bishara:2017pfq . However, at variance with Ref. Bishara:2017pfq , we already find such contributions in the context of the Pauli form factor $F_{1}$ in the $VV$ channel. This is reflected by the corresponding coefficients of $\mathcal{O}_{3}$ and $\mathcal{O}_{5}$. We find similar discrepancies to Ref. Bishara:2017pfq in the NR expansion of the tensor current. Besides the $\mathcal{O}_{3}$ and $\mathcal{O}_{5}$ coefficients, we also disagree in that in our expressions the induced form factors $F_{2,T}$ and $F_{3,T}$ [see Eq. (A) below] combine to the tensor magnetic moments, so that the less well determined individual form factors are not required at this order in the expansion.

Expressed in terms of nucleon matrix elements we have

[TABLE]

and $f_{1}^{V,n}$, $f_{2}^{V,n}$, $\dot{f}_{1}^{V,n}$ are given by the exchange $u\leftrightarrow d$ in $f_{1}^{V,p}$, $f_{2}^{V,p}$, $\dot{f}_{1}^{V,p}$. Finally, $\tilde{f}_{1}^{V,N}$ follows from $f_{1}^{V,N}$ by replacing $C_{q}^{VV}\to C_{q}^{AV}$. For the axial-vector and pseudoscalar matrix elements $g_{A}^{N}$, $g_{P}^{N}$, and $h_{5}^{N}$ we refer to Ref. Hoferichter:2015ipa , since in the present paper the numerical analysis is restricted to the coherently enhanced contributions.

The scalar couplings $f_{q}^{N}$, scalar radii $\dot{\sigma}$ and $\dot{\sigma}_{s}$, as well as charge radii $\langle r_{E}^{2}\rangle^{N}$, $\langle r_{E,s}^{2}\rangle^{N}$ and magnetic moments $\kappa_{N}$, $\kappa_{N}^{s}$ are discussed in detail in Ref. Hoferichter:2016nvd (see also Ref. Hill:2014yxa for the heavy-quark couplings). For strangeness in the vector channel there is increasing evidence from lattice QCD that these couplings are extremely small, we take $\kappa_{N}^{s}=0.006(4)$ and $\langle r_{E,s}^{2}\rangle^{N}=0.0012(9)\,\text{fm}^{2}$ from Ref. Alexandrou:2018zdf . The scalar couplings to $u$- and $d$-quarks can be reconstructed Crivellin:2013ipa from the pion–nucleon $\sigma$-term and input for the proton–neutron mass difference Gasser:1974wd ; Borsanyi:2014jba ; Gasser:2015dwa ; Brantley:2016our . Here, the tension between phenomenology Hoferichter:2015dsa ; Hoferichter:2015hva ; Hoferichter:2016ocj and lattice QCD Durr:2015dna ; Yang:2015uis ; Abdel-Rehim:2016won ; Bali:2016lvx already discussed in Ref. Hoferichter:2016nvd persists. More recently, the phenomenological determination from data on pionic atoms Gotta:2008zza ; Strauch:2010vu ; Hennebach:2014lsa ; Baru:2010xn ; Baru:2011bw has been confirmed by an independent extraction from low-energy pion–nucleon cross sections RuizdeElvira:2017stg , making the resolution of the tension all the more pressing.

The tensor form factors of the nucleon are defined according to Weinberg:1958ut ; Adler:1975he

[TABLE]

where the tensor charges of the proton $f_{q}^{T,p}=F_{1,T}^{q,p}(0)$ have been recently calculated to high precision in lattice QCD (evaluated at scale $\mu=2\,\text{GeV}$) Gupta:2018lvp :

[TABLE]

and the neutron ones follow from isospin symmetry according to

[TABLE]

The other couplings $\tilde{\kappa}_{q}^{T,p}=F_{2,T}^{q,p}(0)+2F_{3,T}^{q,p}(0)$, related to the tensor magnetic moments $\kappa_{q}^{T,p}=-2\tilde{\kappa}_{q}^{T,p}$, are less well determined. Lattice QCD information Gockeler:2006zu is consistent with an estimate based on analyticity and unitarity of the form factors in analogy to Ref. Cirigliano:2017tqn , combining the pion tensor charge Baum:2011rm with the electromagnetic form factors of the nucleon Hoferichter:2016duk . The resulting values $\tilde{\kappa}_{u}^{T,p}=-1.3(5)$, $\tilde{\kappa}_{d}^{T,p}=-0.7(3)$, $\tilde{\kappa}_{s}^{T,p}=0.00(1)$ Hoferichter:2018zwu indicate that for $u$- and $d$-quarks the induced terms in the tensor decomposition in Eq. (A) are actually dominant. As in the vector case, the strangeness content is very small. In addition, the known pion tensor charge Baum:2011rm allows us to calculate the pion-exchange diagram also in this channel, but, similarly to the vector current, this contribution is suppressed by $(N-Z)/A$ due to its isospin structure.

Finally, the spin-$2$ couplings of the nucleon are given as moments of nucleon PDFs $q(x)$ in complete analogy to Eq. (10):

[TABLE]

subject to the same sum rule

[TABLE]

Numerically, this gives for the proton (again at $\mu=2\,\text{GeV}$) Martin:2009iq ; Buckley:2014ana

[TABLE]

while the neutron couplings follow as in Eq. (43) by exchanging $u$ and $d$.

Appendix B Matching to NREFT for scalar dark matter

For a spin-[math], Standard-Model singlet $\chi$ the analysis is based on the effective Lagrangian

[TABLE]

where the notation follows closely the spin-$1/2$ case in Eq. (3), and for a real scalar we have $C_{q}^{VV}=C_{q}^{TT}=0$. We have not included an axial term of the form $\chi^{\dagger}i\partial^{\mu}\chi\,\bar{q}\gamma_{\mu}\gamma_{5}q$, because such a contribution reduces to a combination of $\mathcal{O}_{7}$ and $\mathcal{O}_{10}$ NREFT operators, and is therefore even further suppressed than the standard SD interaction in the spin-$1/2$ case. Similarly, we have neglected a tensor operator with $\bar{q}\sigma_{\mu\nu}i\gamma_{5}q$.

Most single-nucleon amplitudes come out as in the spin-$1/2$ case, up to an additional factor of $m_{\chi}$ whenever a derivative is required in the effective operator, and the corresponding factors of $\Lambda$. To have the reduced amplitudes in the same conventions as for spin-$1/2$, a factor $m_{\chi}$ needs to be removed which leads to

[TABLE]

where we have kept the same notation for the couplings as in Eq. (40), with the understanding that the Wilson coefficients therein now refer to the operators defined in Eq. (47).

In total, the analog of Eq. (II.1) for spin-[math] becomes

[TABLE]

where now $\zeta=1(2)$ for a complex (real) scalar. At this order in the ChEFT expansion we do not find contributions from $\mathcal{O}_{5,8,11}$.

Appendix C Two-body structure factors

The two-body $\theta^{\mu}_{\mu}$ and spin-$2$ amplitudes given in Eqs. (22) and (23) can be rewritten according to

[TABLE]

where we have used that $\mathbf{q}=-\mathbf{q}_{1}-\mathbf{q}_{2}$ and thus

[TABLE]

In both cases most of the first term excluding the $q^{2}$-dependent piece can be related to $\mathcal{M}_{2,\text{NR}}^{SS}$, leading to the $\mathcal{F}_{\pi}$ contribution in Eq. (IV.2). The $q^{2}$-dependent part can be absorbed into a redefinition of $\mathcal{F}_{\text{b}}$ in Eq. (IV.2) to avoid the introduction of another structure factor. This redefinition does not change the normalization of $\mathcal{F}_{\text{b}}(0)$ given by the nuclear binding energy. The remaining term involving pion propagators can be expressed in terms of

[TABLE]

On the other hand, the relativistic corrections to the one-body $\theta$ term and the spin-$2$ channel read

[TABLE]

In the limit $\mathbf{q}\rightarrow 0$ these amplitudes become proportional to the kinetic-energy operator $T$:

[TABLE]

Summarizing all 1b and 2b contributions we obtain

[TABLE]

In the limit $\mathbf{q}\rightarrow 0$ we have $\mathbf{q}_{1}=-\mathbf{q}_{2}$, so that the momentum transfers in $\mathcal{M}_{NN}^{(1)}$ and $\mathcal{M}_{NN}^{(2)}$ become equal and both can be identified with the pion-exchange part $V_{\text{OPE}}^{(0)}$ of the leading-order chiral $NN$ potential:

[TABLE]

Together with the contact terms in Eqs. (57) and (58) we recover twice the complete leading-order chiral $NN$ potential $V_{NN}^{\text{LO}}$.

In the limit $\mathbf{q}\rightarrow 0$ both $\Delta\mathcal{M}_{1,\text{NR}}^{\theta}$ and $\Delta\mathcal{M}_{1,\text{NR}}^{(2)}$ become proportional to $T$, leading in both cases to the linear combination $T+V_{NN}^{\text{LO}}$, that can be renormalized to the nuclear binding energy at LO order in the chiral expansion. For the $\theta^{\mu}_{\mu}$ current in Eqs. (54) and (22) this conclusion follows directly from Eq. (57). On the other hand, the spin-$2$ one-body contribution in Eq. (55) carries a coupling $f_{N}^{(2)}$ that, in general, may differ from $f_{\pi}^{(2)}$, the coupling of the two-body term in Eq. (23), see Eqs. (40) and (II.2). However, due to the sum rules in Eqs. (11) and (45), when all Wilson coefficients are set equal we have $f_{N}^{(2)}/m_{N}=f_{\pi}^{(2)}/M_{\pi}$, so that these two couplings cancel in the last term of Eq. (58) to ensure renormalizability at this order in ChEFT. Since the comparison of Eqs. (II.2) and (A) shows that the individual spin-$2$ couplings do not differ much (especially when considering an isospin average for the nucleon), we assume that for spin-$2$ the kinetic-energy operator aligns as required. Similarly, in Eq. (23) we also assume that the coefficients of the spin-$2$ $(N^{\dagger}N)^{2}$ contact operators can be determined in the same way, at least up to higher-order effects.

In practice, once the contribution at $q=0$ is renormalized to the nuclear binding energy at LO by adjusting the contact terms $C_{S}$ and $C_{T}$, the contributions from $\Delta\mathcal{M}_{1,\text{NR}}^{\theta}$ and $\Delta\mathcal{M}_{1,\text{NR}}^{(2)}$ to the full structure factors are small. In addition, while in principle the structure factors could differ for finite $q$, due to the different functional form of Eqs. (54) and (55) for $q\neq 0$, we find that the results are practically the same using either expression. In view of the uncertainties from higher orders in the chiral expansion, this shows that the definition of a single new structure factor is sufficient. We choose

[TABLE]

which in the limit $\mathbf{q}\rightarrow 0$ reduces to

[TABLE]

The corresponding response functions in the naive shell model are

[TABLE]

where due to Eq. (61) we have $\bar{\mathcal{F}}_{\text{b}}(0)=-2E_{\text{b}}/M_{\pi}$ with the binding energy of the nucleus $E_{\text{b}}<0$. In Eq. (IV.2) we redefine $\mathcal{F}_{\text{b}}(q^{2})=\bar{\mathcal{F}}_{\text{b}}(q^{2})-\frac{q^{2}}{M_{\pi}^{2}}\mathcal{F}_{\pi}(q^{2})$.

We evaluate the two-body structure factors in Eq. (C) in a non-interacting shell model, as described in detail in Ref. Hoferichter:2016nvd , with the multidimensional integrations performed using the CUBA library Hahn:2004fe . For orbitals in the configuration space of the nuclear shell-model calculations, the occupation numbers are taken from the full diagonalization discussed in Sec. III. The final results are represented by a fit function

[TABLE]

For one-body contributions this is an exact analytic form with known maximal power $m$, and the overall normalization is set by $c_{0}=A$ and $c_{0}=Z-N$ for $\mathcal{F}^{M}_{+}$ and $\mathcal{F}^{M}_{-}$, respectively. The expansion is organized in the variable $u=q^{2}b^{2}/2$, with harmonic-oscillator lengths

[TABLE]

that depend on the mass number $A$ of the nucleus. The functional form in Eq. (63) proves an efficient representation of the two-body structure factors as well, but the optimal maximal power needs to be determined empirically. Table 6 provides the complete list.

Appendix D ChiralEFT4DM: Python Notebook

We strongly encourage current and future analyses of direct detection experiments to use the structure factors for the different interaction channels discussed in this manuscript. For convenience we offer an accompanying Python package in the form of a Jupyter notebook, which can be downloaded from https://theorie.ikp.physik.tu-darmstadt.de/strongint/ChiralEFT4DM.html. The notebook calculates both nuclear structure factors and differential recoil spectra:

[TABLE]

including the general response

[TABLE]

where $\rho$ and $f(v)$ denote the local dark matter density and velocity distribution, respectively. The notebook gives results for all possible coherently enhanced couplings of spin-$1/2$ and spin-[math] WIMPs to one and two nucleons up to third order in ChEFT, as discussed in the present manuscript. All stable isotopes of the most relevant nuclear targets including fluorine, silicon, argon, germanium, and xenon are available. In addition, our package calculates the responses based upon the fundamental couplings at the level of quarks and gluons as incorporated in the respective Wilson coefficients.

The notebook aims to be self-explanatory and easy-to-use even for users new to Python. When downloaded from the website the files are stored in an archive. When unpacked to a common directory the notebook can be loaded. In the first part of the notebook, users can specify a given response and create data sets for both nuclear structure factors and differential recoil spectra. In the second part of the notebook, users can set specific values for the Wilson coefficients that describe the WIMP–quark/gluon couplings. The notebook generates the corresponding nucleon and pion matrix elements. Finally, the package yields the response including all channels that contribute to the choice of Wilson coefficients.

For completeness, we also implemented a routine that calculates the rate corresponding to the standard halo model (SHM) Lewin:1995rx . In this way we clarify conventions and facilitate the use of improved astrophysical input in the future, as is increasingly becoming available with the Gaia mission Brown:2018dum , see, e.g., Refs. Necib:2018iwb ; Eilers:2018 ; Evans:2018bqy . The SHM is defined by $\rho=0.3\,\text{GeV}/\text{cm}^{3}$ and

[TABLE]

where $z=v_{\text{esc}}/v_{0}$, $v_{\text{esc}}=544\,\text{km}/\text{s}$, $v_{0}=220\,\text{km}/\text{s}$, and the Earth’s velocity $v_{E}=232\,\text{km}/\text{s}$ drops out in the normalization. For the operators considered in the notebook one needs

[TABLE]

where $x_{\text{min}}=v_{\text{min}}/v_{0}$, $\eta=v_{E}/v_{0}$, and

[TABLE]