The computation of massively separated flows is a challenging problem of particular in- terest in industrial applications. For the purpose of properly simulating these complex flows on not too heavy unstructured meshes as usually employed in industry, appropriate numerical and turbulent models must be used. In the present work, the computation of the flow past a circular cylinder at different Reynolds numbers is chosen as benchmark. The spatial discretization is based on a mixed finite element/finite volume formulation on unstructured grids. The numerical dissipation of the upwind scheme is made of sixth-order space derivatives in order to limit as far as possible the interactions between numerical and subgrid scale (SGS) dissipation, which could deteriorate the accuracy of the results [4]. A variational multi-scale large-eddy simulation (VMS-LES) with dynamic SGS models and a RANS/VMS-LES model are evaluated on the proposed benchmark for subcritical and supercritical flow regimes respectively (see Fig. 1 and Tab. 1). In the VMS-LES used in this work, the separation between the large and the small resolved scales is obtained through a variational projection operator based on spatial average on agglomerated cells [1]. The dynamic procedure allows the adaptation of the constant of the SGS model to the spatial and temporal variation of the flow characteristics, while the VMS formulation restricts the SGS model effects to the smallest resolved scales. The dynamic versions of the Smagorinsky and of the WALE SGS models are considered herein. The non-dynamic counterparts of these SGS models are also used in order to evaluate the impact of dy- namic SGS modeling in the considered VMS-LES approach for the simulation of massively separated flows. However, the Reynolds number range useful for LES-like simulation is limited as LES grid needs to be sufficiently fine to resolve a significant part of the turbulence scales. With the aim of simulating high Reynolds number flows, it is considered in the present work a hybridization strategy using a blending parameter, such that a VMS- LES simulation is obtained where the grid resolution is fine enough to resolve a significant part of the turbulence fluctuations [2], while a RANS model is acting in the regions of coarse grid resolution, as, for instance, near the body surface