Determining the leading-order contact term in neutrinoless double $\boldsymbol{\beta}$ decay

We present a method to determine the leading-order (LO) contact term contributing to the $nn \to pp e^-e^-$ amplitude through the exchange of light Majorana neutrinos. Our approach is based on the representation of the amplitude as the momentum integral of a known kernel (proportional to the neutrino propagator) times the generalized forward Compton scattering amplitude $n(p_1) n(p_2) W^+ (k) \to p(p_1^\prime) p(p_2^\prime) W^- (k)$, in analogy to the Cottingham formula for the electromagnetic contribution to hadron masses. We construct model-independent representations of the integrand in the low- and high-momentum regions, through chiral EFT and the operator product expansion, respectively. We then construct a model for the full amplitude by interpolating between these two regions, using appropriate nucleon factors for the weak currents and information on nucleon-nucleon ($N\! N$) scattering in the $^1S_0$ channel away from threshold. By matching the amplitude obtained in this way to the LO chiral EFT amplitude we obtain the relevant LO contact term and discuss various sources of uncertainty. We validate the approach by computing the analog $I = 2$ $N\! N$ contact term and by reproducing, within uncertainties, the charge-independence-breaking contribution to the $^1S_0$ $N\! N$ scattering lengths. While our analysis is performed in the $\overline{\rm MS}$ scheme, we express our final result in terms of the scheme-independent renormalized amplitude ${\cal A}_\nu(|{\bf p}|,|{\bf p}^\prime|)$ at a set of kinematic points near threshold. We illustrate for two cutoff schemes how, using our synthetic data for ${\cal A}_\nu$, one can determine the contact-term contribution in any regularization scheme, in particular the ones employed in nuclear-structure calculations for isotopes of experimental interest.