A Galerkin finite element method is considered to approximate the incompressible Navier-Stokes equations together with iterative methods to solve a resulting system of algebraic equations. This system couples velocity and pressure unknowns, thus requiring a special technique for handling. We consider the Navier-Stokes equations in velocity--kinematic pressure variables as well as in velocity--Bernoulli pressure variables. The latter leads to the rotation form of nonlinear terms. This form of the equations plays an important role in our studies. A consistent stabilization method is considered from a new view point. Theory and numerical results in the paper address both the accuracy of the discrete solutions and the effectiveness of solvers and a mutual interplay between these issues when particular stabilization techniques are applied.