  We prove from the modified gravity (MG) galaxy/halo mock catalogs that the
shape of the pairwise kinematic Sunyaev-Zeldovich (kSZ) power spectrum $P_{\rm
kSZ}$ has constraining power on discriminating different gravity theories on
cosmological scales. By varying the effective optical depth $\tau_{\rm T}$ as a
free parameter, we verify that the $\tau_{\rm T}$-$f$ (the linear growth rate)
degeneracy in the linear theory of $P_{\rm kSZ}$ is broken down by the
non-linear structure growth and the scale-dependence of $f$ in some MG
theories. Equivalently speaking, the shape of $P_{\rm kSZ}$ alone could be used
to constrain the MG theories on cosmological scales. With a good knowledge of
galaxy density biases, we verify that, a combination of the next generation of
CMB and galaxy spectroscopic redshift surveys, e.g. BOSS+CMB-S4 or DESI+CMB-S4,
could potentially discriminate $f(R)$ models from the general relativity at
$\sim5\sigma$ level using the shape of the galaxy pairwise kSZ dipole $P_{{\rm
kSZ},\ell=1}$ alone, when $f_{R0}=10^{-4}$.
