Axions can also couple to the electric dipole moment of nucleons through the operator $-i{g}_{d}a\overline{N}{\sigma }_{\mu \nu }{\gamma }^{5}N{F}^{\mu \nu }/2$. However, this coupling is stringently constrained by astrophysical observations, and we estimate that the sensitivities of the configurations discussed in this work are insufficient to surpass these constraints.

This gives $\mathit{\sigma }(t)=e^{-i{\omega }_{0}t}(\widehat{\mathbf{x}}\mp i\widehat{\mathbf{y}})$, which corresponds to the spins rotating in the clockwise (or counterclockwise) direction around $\widehat{\mathbf{z}}$ when viewed from above if $-$ (or $+$) is taken on, consistent with the dynamics of spins described by the Bloch equation for $\gamma >0$ (or $<0$). Note that in this complex representation, $\mathrm{Re}[\mathit{\sigma }(t)]$ is a unit vector.

For axion dark matter, the field amplitude is $a=\sqrt{2{\rho }_{\mathrm{DM}}}/m_a\sim {10}^3\mathrm{eV}(m_a/1\mathrm{\mu }\mathrm{eV}{)}^{-1}$ and the field gradient is $|\nabla a|\sim m_a{v}_0{a}_0\sim {10}^{-6}{\mathrm{eV}}^2$. Additionally, $|\nabla a|/a\sim m_a{v}_0\sim {10}^{-9}\mathrm{eV}(m_a/1\mathrm{\mu }\mathrm{eV})$. In these numerical expressions, we assumed an expected local density $\simeq 0.3\mathrm{GeV}/\mathrm{c}{\mathrm{m}}^3$ and a characteristic velocity of $v_0\sim {10}^{-3}$.

The $Q$-factor of NMR for a polarized sample with $T_2=100\mathrm{s}$ and Larmor frequency ${\omega }_0\sim \mathrm{GHz}$ is $Q\equiv {\omega }_0{T}_2/2\sim {10}^9$, resulting in a resonance width $\mathrm{\Delta }\omega /{\omega }_0=Q^{-1}\sim {10}^{-9}$. The longitudinal field ${\mathbf{B}}_0$ for the source and the sensor thus needs to be matched with a relative precision of ${10}^{-9}$ to satisfy the resonance condition for the sensor medium. Current superconducting magnets with $\mathcal{O}(10)\mathrm{T}$ magnetic field can achieve ${10}^{-9}$ uniformity within centimeter-scale volume with ${10}^{-8}/\mathrm{hour}$ shift [44], and ${10}^{-6}$ uniformity within decimeter-scale volume with ${10}^{-9}/\mathrm{hour}$ shift [45]. Persistent efforts are being made to develop more GHz-class NMR magnets with field inhomogeneity below ${10}^{-9}$ to resolve the narrow spectral lines in structure biology [44]. Operating at a lower ${\omega }_0$ can relax the requirement on field homogeneity. As we will discuss in the next section, reducing ${\omega }_0$ does not degrade our projected sensitivity.

The ${\mathbf{B}}_0$ field applied to the spin source can magnetize the shielding material, which induces a residual field, altering the total longitudinal field. This effect is negligible as long as the ratio of the residual field to $B_0$ remains below $Q^{-1}$. Otherwise, it can be solved by introducing an identical magnetic shield for the spin sensor or tuning the longitudinal field with a relative precision of $\mathcal{O}(Q^{-1})$.

Prolonged spin precession has be achieved for $\mathcal{O}({10}^4)\mathrm{s}$ in ${}^{129}\mathrm{Xe}$ and ${}^3\mathrm{He}$ gases using spin masers, which utilize optical pumping along the $z$ direction and a rotating transverse magnetic field with a 90 ° phase delay (compared to the spin precession) generated by a feedback system [48, 49, 50, 51]. By appropriately choosing the pumping and feedback parameters, the stabilized transverse magnetization can achieve $\mathcal{O}(1)$ the initial magnetization $M_0$ [52]. We hope in the future this technique can be further applied to ${}^{129}\mathrm{Xe}$ or ${}^3\mathrm{He}$ liquids.

For the $R=10\mathrm{cm}$ setup, the axion field in the spin sensor slightly deviates from the near-field result for small masses.

Note that the actual magnetic field that fluxes through a pickup loop near the end of the sensor medium is $\mathcal{O}(1)$ smaller than $M_x$ (see Supplemental Material for details). Since $M_x\propto g_{aN}^2$, our sensitivities obtained by treating them as equal are slightly overestimated by an $\mathcal{O}(1)$ factor.

In the $m_a\gg 1/R$ regime, the magnetization within the detection sample decreases dramatically and the signal is only contributed by a small portion of the sample around the near end. In Fig. 1 (right), we always take $V$ in $P_k^{\mathrm{SP}}$ as the volume of the entire detection sample such that the spin projection noise is underestimated in this regime. The actual sensitivities are thus decaying more rapidly when $m_a\gg 1/R$ than in Fig. 1 (right).

This still holds true with the Bayesian framework we adopt, even though our signal is deterministic rather than stochastic. We provide a detailed discussion of noise sources and the likelihood analysis in Supplemental Material.