The answer is no if $k$ is independent on $\lambda$. Take a $2$-dimensional torus with a flat metric, which is a very nice manifold. Eigenfunctions can be found explicitly, and they are analytic. There can be some points where many nodal lines intersect. Near these points, eigenfunctions will look like $u(z)=\Re z^n$ where $n$ can be large. So for every $k$ we can we can find $n$ such that the integral of $1/|u|^k$ diverges. But it looks plausible that the integral will converge with some $k$ depending on $\lambda$ (and on the manifold).
For example, if $\lambda$ is the lowest eigenvalue, then the nodal set is empty, and $1/|u|^k$ is bounded thus summable for every $k$. If $\lambda$ is second lowest, then the nodal set is probably consists of non-singular curves, and $1/|u|^{k}$ is summable for $k<1$. How complex can be a singularity of the nodal set for large eigenvalues I do not know, but perhaps for a torus or a sphere this can be determined.
EDIT. H. Donnelly and Ch. Fefferman proved that an a zero of an eigenfunction can have multiplicity at most $C\sqrt{\lambda}$, where $\lambda$ is the eigenvalue, and $C$ depends on the manifold. Up to the value of $C$, this is best possible: the example is spherical harmonics. This solves the question (up to this constant $C$, with $k$ dependent on $\lambda$.