  We perform a detailed study of the relaxation towards equilibrium in the
Hamiltonian Mean-Field (HMF) model, a prototype for long-range interactions in
$N$-particle dynamics. In particular, we point out the role played by the
infinity of stationary states of the associated $N ~$ Vlasov dynamics. In this
context, we derive a new general criterion for the stability of any spatially
homogeneous distribution, and compare its analytical predictions with numerical
simulations of the Hamiltonian, finite $N$, dynamics. We then propose and
verify numerically a scenario for the relaxation process, relying on the Vlasov
equation. When starting from a non stationary or a Vlasov unstable stationary
initial state, the system shows initially a rapid convergence towards a stable
stationary state of the Vlasov equation via non stationary states: we
characterize numerically this dynamical instability in the finite $N$ system by
introducing appropriate indicators. This first step of the evolution towards
Boltzmann-Gibbs equilibrium is followed by a slow quasi-stationary process,
that proceeds through different stable stationary states of the Vlasov
equation. If the finite $N$ system is initialized in a Vlasov stable homogenous
state, it remains trapped in a quasi-stationary state for times that increase
with the nontrivial power law $N^{1.7}$. Single particle momentum distributions
in such a quasi-stationary regime do not have power-law tails, and hence cannot
be fitted by the $q$-exponential distributions derived from Tsallis statistics.
