  Under some assumptions on the hierarchy of relevant energy scales, we compute
the nonrelativistic QCD (NRQCD) long-distance matrix elements (LDMEs) for
inclusive production of $J/\psi$, $\psi(2S)$, and $\Upsilon$ states based on
the potential NRQCD (pNRQCD) effective field theory. Based on the pNRQCD
formalism, we obtain expressions for the LDMEs in terms of the quarkonium
wavefunctions at the origin and universal gluonic correlators, which do not
depend on the heavy quark flavor or the radial excitation. This greatly reduces
the number of nonperturbative unknowns and substantially enhances the
predictive power of the nonrelativistic effective field theory formalism. We
obtain improved determinations of the LDMEs for $J/\psi$, $\psi(2S)$, and
$\Upsilon$ states thanks to the universality of the gluonic correlators, and
obtain phenomenological results for cross sections and polarizations at large
transverse momentum that agree well with measurements at the LHC.
