We present a new method for modelling the kinematics of galaxies from interferometric observations by performing the optimization of the kinematic model parameters directly in visibility space instead of the conventional approach of fitting velocity fields produced with the CLEAN algorithm in real-space. We demonstrate our method on Atacama Large Millimeter/submillimeter Array (ALMA) observations of $^{12}$CO (2-1), (3-2), or (4-3) emission lines from an initial sample of 30 massive 850 $\mu$m-selected dusty star-forming galaxies with far-infrared luminosities $\gtrsim$$\, 10^{12}$ L$_{\odot }$ in the redshift range $z \sim$ 1.2-4.7. Using the results from our modelling analysis for the 12 of the 20 sources with the highest signal-to-noise emission lines that show disc-like kinematics, we conclude the following: (i) our sample prefers a CO-to-$H_2$ conversion factor, of $\alpha _{\rm CO} = 0.74 \pm 0.37$; (ii) these far-infrared luminous galaxies follow a similar Tully-Fisher relation between the circular velocity, $V_{\rm circ}$, and baryonic mass, $M_{\rm b}$, as less strongly star-forming samples at high redshift, but extend this relation to much higher masses - showing that these are some of the most massive disc-like galaxies in the Universe; (iii) finally, we demonstrate support for an evolutionary link between massive high-redshift dusty star-forming galaxies and the formation of local early-type galaxies using the both the distributions of the baryonic and kinematic masses of these two populations on the $M_{\rm b}$ - $\sigma$ plane and their relative space densities.