  In this paper, we study the Sobolev regularity of solutions to nonlinear
second order elliptic equations with super-linear first-order terms on
Riemannian manifolds, complemented with Neumann boundary conditions, when the
source term of the equation belongs to a Lebesgue space, under various
integrability regimes. Our method is based on an integral refinement of the
Bochner's identity, and leads to "semilinear Calder\'on-Zygmund" type results.
Applications to the problem of smoothness of solutions to Mean Field Games
systems with Neumann boundary conditions posed on convex domains of the
Euclidean space will also be discussed.
