  Massive neutrinos suppress the growth of structure on small scales and leave
an imprint on large-scale structure that can be measured to constrain their
total mass, $M_\nu$. With standard analyses of two-point clustering statistics,
$M_\nu$ constraints are severely limited by parameter degeneracies. Hahn et
al.(2020) demonstrated that the bispectrum, the next higher-order statistic,
can break these degeneracies and dramatically improve constraints on $M_\nu$
and other cosmological parameters. In this paper, we present the constraining
power of the redshift-space galaxy bispectrum, $B_0^g$. We construct the Molino
suite of 75,000 mock galaxy catalogs from the Quijote $N$-body simulations
using the halo occupation distribution (HOD) model, which provides a galaxy
bias framework well-suited for simulation-based approaches. Using these mocks,
we present Fisher matrix forecasts for $\{\Omega_m,\Omega_b,h,n_s,\sigma_8,
M_\nu\}$ and quantify, for the first time, the total information content of
$B_0^g$ down to nonlinear scales. For $k_{\rm max}=0.5h/Mpc$, $B_0^g$ improves
constraints on $\Omega_m,\Omega_b,h,n_s,\sigma_8$, and $M_\nu$ by 2.8, 3.1,
3.8, 4.2, 4.2, and $4.6\times$ over the power spectrum, after marginalizing
over HOD parameters. Even with priors from $Planck$, $B_0^g$ improves all of
the cosmological constraints by $\gtrsim2\times$. In fact, for $P_\ell^g$ and
$B_0^g$ out to $k_{\rm max}=0.5h/Mpc$ with $Planck$ priors, we achieve a
$1\sigma$ $M_\nu$ constraint of 0.048 eV, which is tighter than the current
best cosmological constraint. While effects such as survey geometry and
assembly bias will have an impact, these constraints are derived for
$(1h^{-1}{\rm Gpc})^3$, a substantially smaller volume than upcoming surveys.
Therefore, we conclude that the galaxy bispectrum will significantly improve
cosmological constraints for upcoming galaxy surveys -- especially for $M_\nu$.
