Ultra-luminous X-ray sources (ULXs) are X-ray sources located off-galactic centre and have luminosities exceeding the Eddington limit of a stellar-mass black hole ($L_X>10^{39}\;{\rm erg\,s}^{-1}$). With the discovery of X-ray pulsations in some of these objects (e.g. M82 X-2), we know a part of the ULX population has a neutron star (NS) accretor. We present systematic modelling of intermediate-mass X-ray binaries (IMXBs; donor-star mass range $2.0$ to $8.0\;$M$_{\odot}$) to explain the formation of this sub-population of ULXs. Using MESA, we explore the allowed initial parameter space of binary systems consisting of a neutron star and an intermediate-mass donor star, that could explain the observed properties of ULXs. Our simulations take into account beaming effects and also include stellar rotation, tides, general angular momentum losses, and a detailed and self-consistent calculation of the mass-transfer rate. Exploring the initial parameters that lead to the formation of neutron-star ULXs, we study the conditions that lead to dynamical stability of these systems, which depends strongly on the response of the donor star to mass loss. Using two values for the initial NS mass ($1.3\;$M$_{\odot}$ and $2.0\;$M$_{\odot}$) we find that IMXBs can produce NS-ULXs with typical time-averaged isotropic-equivalent X-ray luminosities between $10^{39}$--$10^{41}\;{\rm erg\,s}^{-1}$ on a timescale of up to $\sim\!1.0\;$Myr for the lower luminosities, assuming geometrical beaming. Finally, we estimate their likelihood of detection, the types of white-dwarf remnants left behind by the donors, as well as the total amount of mass accreted by the neutron stars. We also compare our results to the observed pulsating ULXs and infer their initial parameters. Our results suggest that a large subset of the observed pulsating ULX population can be explained by IMXBs in a super-Eddington mass-transfer phase.