Fragmentation is the dominant mechanism for hadron production with high transverse momentum. For spin-triplet S-wave heavy quarkonium production, contribution of gluon fragmenting to color-singlet channel has been numerically calculated since 1993. However, there is still no analytic expression available up to now because of its complexity. In this paper, we calculate both polarization-summed and polarized fragmentation functions of gluon fragmenting to a heavy quark-antiquark pair with quantum number $^3S_1^{[1]}$. Our calculations are performed in two different frameworks. One is the widely used nonrelativistic QCD factorization, and the other is the newly proposed soft gluon factorization. In either case, we calculate at both leading order and next-to-leading order in velocity expansion. All of our final results are presented in terms of compact analytic expressions.