As a consequence of the Helmholtz–Hodge theorem, any divergence-free vector field can be decomposed in two L 2-orthogonal, solenoidal vector fields expressed in terms of projections of the velocity and vorticity fields, on an arbitrary direction in space. Based on this type of decomposition and the choice of the wall-normal direction, an efficient spectral code is developed for incompressible flows developing between two parallel walls. The method relies on a weak formulation of the Navier–Stokes equations in the two corresponding divergence-free subspaces. The approximation is based on Fourier expansions in two directions and on the Chebyshev basis proposed by Moser et al. in the third direction in order to satisfy the wall boundary conditions. The method accuracy is validated for the plane Poiseuille linear stability problem and compared with the case of a spectral collocation method. Simulations of by-pass transition in boundary layers developing between two parallel walls are then presented. Since, by construction, the two orthogonal vector fields of the decomposition are associated respectively to the Orr–Sommerfeld and the Squire modes of the linear stability theory, the method makes it possible to evaluate kinetic energy transfers due to the coupling between these two scalar modes and their interactions with the base flow. The decomposition is also used to describe the structure of finite-length streaks in the earlier stages of transition.