We study the generation of vorticity and velocity dispersion by orbit crossing using cosmological numerical simulations, and calculate the backreaction of these effects on the evolution of large-scale density and velocity divergence power spectra. We use Delaunay tessellations to define the velocity field, showing that the power spectra of velocity divergence and vorticity measured in this way are unbiased and have better noise properties than for standard interpolation methods that deal with mass weighted velocities. We show that high resolution simulations are required to recover the correct large-scale vorticity power spectrum, while poor resolution can spuriously amplify its amplitude by more than one order of magnitude. We measure the scalar and vector modes of the stress tensor induced by orbit crossing using an adaptive technique, showing that its vector modes lead, when input into the vorticity evolution equation, to the same vorticity power spectrum obtained from the Delaunay method. We incorporate orbit crossing corrections to the evolution of large scale density and velocity fields in perturbation theory by using the measured stress tensor modes. We find that at large scales (k ≃ 0.1 h Mpc −1 ) vector modes have very little effect in the density power spectrum, while scalar modes (velocity dispersion) can induce percent level corrections at z = 0, particularly in the velocity divergence power spectrum. In addition, we show that the velocity power spectrum is smaller than predicted by linear theory until well into the nonlinear regime, with little contribution from virial velocities.