Summary: We consider the numerical solution of quasilinear elliptic Neumann problems \[ \begin{cases} T(u) &\equiv -\text{div}(f(x,\nabla u)\nabla u)= g(x),\\ {\partial u\over\partial\nu}|_{\partial\Omega} &= 0.\end{cases} \] The basic difficulty is the non-injectivity of the operator, which can be overcome by suitable factorization. We extend the gradient-finite element method, introduced earlier by the authors for Dirichlet problems, to the Neumann problem. The algorithm is constructed and its convergence is proved.