  We consider a generic dispersive massive gravity theory and numerically study
its resulting modified energy and strain spectra of gravitational waves (GWs)
sourced by (i) fully developed turbulence during the electroweak phase
transition (EWPT); and (ii) forced hydromagnetic turbulence during the QCD
phase transition (QCDPT). The GW spectra are then computed in both the spatial
and temporal Fourier domains. We find, from the spatial spectra, that the slope
modifications are independent of the eddy size at QCDPT, and, from the temporal
spectra, that the modifications are pronounced in the nHz-10nHz range -- the
sensitivity range of the North American Nanohertz Observatory for Gravitational
Waves (NANOGrav) -- for a graviton mass $m_{\rm g}$ in the range
$2.2\times10^{-23}{\rm eV}/c^2\lesssim m_{\rm g}\lesssim7.4\times10^{-22}{\rm
eV}/c^2$.
