Three-gluon decays of radially excited quarkonia $\psi (2S)$ and $\mathrm{{\rm Y}}(2S)$ with both relativistic and QCD radiative corrections

For radially excited heavy quarkonia $V=\psi (2S)$ and $\mathrm{{\rm Y}}(2S)$, the nodal structure of their wave function renders the three-gluon decay $V\to ggg$ acutely sensitive to relativistic corrections, presenting a longstanding challenge for reliable theoretical predictions. We perform a comprehensive analysis of these decays within the Bethe-Salpeter formalism, constructing analytic harmonic oscillator wave functions that explicitly incorporate the $2S$ nodal structure. Model-independent relations among polarized decay widths are derived from helicity-flip and phase-space symmetries. To obtain physically consistent results beyond $q^2$-order relativistic corrections, we introduce a concise phenomenological treatment that effectively incorporates partial higher-order contributions while preserving the correct low-momentum limit. Taking into account both relativistic and QCD radiative corrections, we compute $\mathrm{\Gamma }(V\to ggg)$, $\mathrm{\Gamma }(V\to e^{+}{e}^{-})$, and their ratio $R_V=\mathrm{\Gamma }(V\to ggg)\mathrm{/}\mathrm{\Gamma }(V\to e^{+}{e}^{-})$, which allows a clean extraction of the harmonic oscillator parameter ${\beta }_V$, with the extracted values lying at the lower end of the ranges from typical phenomenological models. Our predictions for the branching ratios $\mathcal{B}(V\to ggg)$ and $\mathcal{B}(V\to e^{+}{e}^{-})$ are in excellent agreement with experimental data. A striking feature is the distinct convergence behavior of the relativistic expansion: the leptonic width converges rapidly already at $q^2$ order, whereas the gluonic width converges far more slowly due to destructive interference induced by the nodal structure in the multi-gluon convolution. These results offer new insights into the dynamics of excited quarkonia and provide valuable constraints on non-perturbative descriptions of heavy-quark bound states.