  Recently, the first successful attempt at computing stellar models in two
dimensions has been presented with models that include the centrifugal
deformation and self-consistently compute the velocity field. This paper aims
at studying the rotational evolution of 2D models of stars rotating at a
significant fraction of their critical angular velocity. From the predictions
of these models, we aim to improve our understanding of the formation of single
Be stars. Using the ESTER code that solves the stellar structure of a rotating
star in two dimensions with time evolution, we have computed evolution tracks
of stars between 4 and 10Msun for initial rotation rates ranging between 60 and
90% the critical rotation rate. A minimum initial rotation rate at the start of
the main sequence is needed to spin up the star to critical rotation within its
main sequence lifetime. This threshold depends on the stellar mass, and
increases with increasing mass. The models do not predict any stars above 8Msun
to reach (near) critical rotation during the main sequence. Furthermore, we
find the minimum threshold of initial angular velocity is lower for SMC
metallicity compared to Galactic metallicity, which is in agreement with the
increased fraction in the number of observed Be stars in lower metallicity
environments. Self-consistent 2D stellar evolution provide more insight into
the rotational evolution of intermediate-mass stars, and our predictions are
consistent with observations of velocity distributions and fraction of Be stars
amongst B-type stars. We find that stars with a mass above 8Msun do not
increase their fraction of critical rotation during the main sequence. Since a
fraction of stars above 8Msun have been observed to display the Be phenomenon,
other processes or formation channels must be at play, or critical rotation is
not required for the Be phenomenon above this mass.
