The aim of this paper is to study the properties of approximations to nonlinear terms
of the 2D incompressible Navier-Stokes equations in the stream function formulation (time-dependent
biharmonic equation). The nonlinear convective terms are numerically solved by
using the method with internal iterations, compared to the ones which are solved by using
explicit and implicit schemes (operator splitting scheme Christov and Marinova; (2001)). Using schemes and algorithms,
the steady 2D incompressible flow in a lid-driven cavity is solved up to Reynolds number Re
=5000 with second-order spatial accuracy. The schemes are thoroughly validated on
grids with different resolutions. The result of numerical experiments shows that the finite difference
scheme with internal iterations on nonlinearity is more efficient for the high Reynolds number.