  Within the nonrelativistic QCD (NRQCD) factorization framework, we compute
the ${\mathcal O}(\alpha_s v^2)$ corrections to the hadronic decay rate of
vector quarkonia, exemplified by $J/\psi$ and $\Upsilon$. Setting both the
renormalization and NRQCD factorization scales to be $m_Q$, we obtain
$\Gamma(J/\psi\to {\rm LH})= 0.0716\frac{\alpha_s^3}{m_c^2} \langle
\mathcal{O}_1({}^3S_1)\rangle_{J/\psi}
[1-1.19\alpha_s+(-5.32+3.03\alpha_s)\langle v^2\rangle_{J/\psi}]$ and
$\Gamma(\Upsilon\to {\rm LH})=
0.0716\frac{\alpha_s^3}{m_b^2}\langle\mathcal{O}_1({}^3S_1)\rangle_{\Upsilon}[1-1.56\alpha_s+(-5.32+4.61\alpha_s)\langle
v^2\rangle_{\Upsilon}]$. We confirm the previous calculation of
$\mathcal{O}(\alpha_s)$ corrections on a diagram-by-diagram basis, with the
accuracy significantly improved. For $J/\psi$ hadronic decay, we find that the
${\mathcal O}(\alpha_sv^2)$ corrections are moderate and positive, nevertheless
unable to counterbalance the huge negative corrections. On the other hand, the
effect of ${\mathcal O}(\alpha_sv^2)$ corrections for $\Upsilon(nS)$ is
sensitive to the $\mathcal{O}(v^2)$ NRQCD matrix elements. With the appropriate
choice of the NRQCD matrix elements, our theoretical predictions for the decay
rates may be consistent with the experimental data for $\Upsilon(1S,2S)\to {\rm
LH}$. As a byproduct, we also present the theoretical predictions for the
branching ratio of $J/\psi(\Upsilon)\to 3\gamma$ accurate up to
$\mathcal{O}(\alpha_s v^2)$.
