Abstract. A family of virtual element methods for the two-dimensional Navier--Stokes equations is proposed and analyzed. The schemes provide a discrete velocity field which is pointwise divergence-free. A rigorous error analysis is developed, showing that the methods are stable and optimally convergent. Several numerical tests are presented, confirming the theoretical predictions. A comparison with some mixed finite elements is also performed.Key words. virtual element method, polygonal meshes, Navier--Stokes equations AMS subject classifications. 65N30, 76D05 DOI. 10.1137/17M11328111. Introduction. The virtual element method (VEM), introduced in [11,12], is a recent paradigm for the approximation of partial differential equation problems that shares the same variational background as the finite element methods. The original motivation of VEM is the need to construct an accurate conforming Galerkin scheme with the capability to deal with highly general polygonal/polyhedral meshes, including``hanging vertexes"" and nonconvex shapes. Among the Galerkin schemes, VEM is peculiar in that the discrete spaces consist of functions which are not known pointwise, but about which a limited set of information is available. This limited information is sufficient to construct the stiffness matrix and the right-hand side.The VEM has been developed for many problems; see, for example, [23,1,10,43,46,9,19,17,18,6,42,54,51,37,45]. More specifically, with regard to the Stokes problem, virtual elements have been developed in [3,28,15,25,26,52]. Moreover, VEM is also attracting growing interest for continuum mechanics problems within the engineering community. We cite here the recent works [36,13,4,53,30,2,35] and [8,24,5], for instance. Finally, some examples of other numerical methods for the Stokes or Navier--Stokes equations that can handle polytopal meshes are [33,47,32].In this paper, we initiate the development of the VEM for the Navier--Stokes equations. We limit the study to two-dimensional domains and to diffusion dominated cases. Although this is the simplest situation, it is nonetheless challenging and shows the satisfactory performance of the method; considering more complex cases will be a further step in future work. The presented scheme may be considered as a natural evolution of our recent divergence-free approach developed in [15] for the Stokes problem. However, the nonlinear convective term in the Navier--Stokes equations leads to the introduction of suitable projectors. These, in turn, suggest making use of an enhanced discrete velocity space [52], which is an improvement with respect