Abstract We present new algorithms for the numerical approximation of eigenvalues and invariant subspaces of matrices with cheap action (for example, large but sparse). The methods work with inexact solutions of generalized algebraic Riccati equations. The simpler ones are variants of Subspace Iteration and Block Rayleigh Quotient Iteration in which updates orthogonal to current approximations are computed. Subspace acceleration leads to more sophisticated algorithms. Starting with a Block Jacobi Davidson algorithm, we move towards an algorithm that incorporates Galerkin projection of the non-linear Riccati equation directly, extending ideas of Hu and Reichel in the context of Sylvester equations. Numerical experiments show that this leads to very a competitive algorithm, which we will call the Riccati method, after J.F. Riccati (1676–1754).