Observations of pulsar glitches have the potential to provide constraints on the dynamics of the high density interior of neutron stars. However, to do so, realistic glitch models must be constructed and compared to the data. We take a step towards this goal by testing non-linear models for the mutual friction force, which is responsible for the exchange of angular momentum between the neutron superfluid and the observable normal component in a glitch. In particular, we consider a non-linear dependence of the drag force on the relative velocity between superfluid vortices and the normal component, in which the contributions of both kelvon and phonon excitations are included. This non-linear model produces qualitatively new features, and is able to reproduce the observed bimodal distribution of glitch sizes in the pulsar population. The model also suggests that the differences in size distributions in individual pulsars may be due to the glitches being triggered in regions with different pinning strengths, as stronger pinning leads to higher vortex velocities and a qualitatively different mutual friction coupling with respect to the weak pinning case. Glitches in pulsars that appear to glitch quasi-periodically with similar sizes may thus be due to the same mechanisms as smaller events in pulsars that have no preferred glitch size, but simply originate in stronger pinning regions, possibly in the core of the star.