We investigate the old problem of the fast relaxation of collisionless $N$-body systems which are collapsing or perturbed, emphasizing the importance of (non-collisional) discreteness effects. We integrate orbit ensembles in fixed external potentials, estimating the entropy of the ensemble to analyze the time evolution of the distribution function. We show that these estimates capture the correct physical behavior expected from the 2nd Law of Thermodynamics, without any spurious entropy production. For self-consistent (i.e. stationary) samples, the entropy is conserved, while for non-self-consistent samples, it increases within a few dynamical times up to a maximum where it stabilizes (even in integrable potentials). Our results shed light on the main ingredients for this fast collisionless relaxation. The fundamental ingredient is the discreteness (finite $N$) of gravitational systems in any potential. Additionally, in non-integrable potentials, the presence of chaotic orbits accelerates the entropy production. Contrary to the traditional violent relaxation scenario, our results indicate that a time-dependent potential is not necessary for this fast relaxation. For the first time, in connection with the Nyquist-Shannon theorem we determine the $N$-dependence of this discreteness-driven relaxation, deriving a typical relaxation time $T/\tau_{cr}\approx 0.1 N^{1/6}$, with slightly weaker $N$-dependencies for non-integrable potentials with substantial fractions of chaotic orbits. This timescale is much smaller than the collisional timescale even for small-$N$ systems such as open clusters and represents an upper limit for real collisionless $N$-body systems. Additionally, our results reinforce the conclusion draw in Beraldo e Silva et al. (2017) that the Vlasov equation does not provide an adequate kinetic description of the fast relaxation of collapsing collisionless $N$-body systems.