Large-eddy simulations of atmospheric boundary layer flows have historically relied on the open-channel flow setup, using mixed pseudospectral and finite differences methods. Single-domain spectral solvers are naturally suited for the simulation of homogeneous turbulence in simply-connected domains. Their quality, however, is expected to degrade significantly in the presence of complex geometries or non-periodic flow setups. On the contrary, general-purpose finite volumes solvers are expected to provide a more consistent performance across a range of flow configurations. In view of future applications involving flow over complex terrain, the quality and reliability of a colocated, unstructured, finite-volume solver (OpenFOAM framework) is here analyzed for the simulation of a pressure-driven atmospheric boundary layer flow. A sensitivity analysis on grid resolution, aspect ratio, and CFL condition is performed, and different numerical schemes are also considered. First and second order statistics, as well as velocity spectra and two-point velocity correlations, are compared to predictions from a “battle-tested” pseudo-spectral solver. The solution is found to be particularly sensitive to the grid aspect ratio and to the chosen numerical scheme. First and second order statistics obtained using a non-dissipative setup compare well between the solvers, with the finite volume one featuring an overdissipative behavior that leads to enhanced sub-grid scale stress contributions. When considering velocity spectra, the finite volume solver features a rapid decay of energy density within the inertial subrange, irrespectively of the discretization scheme that is adopted. In addition, a spurious pile up of energy density at high wave numbers is observed across all of the considered cases. The cause of this behavior will be discussed and mitigation strategies proposed.