Let \(\Omega \subset {\mathbb {R}}^N\) ( \(N \ge 3\) ) be a \(C^2\) bounded domain and \(\Sigma \subset \Omega \) is a \(C^2\) compact boundaryless submanifold in \({\mathbb {R}}^N\) of dimension k, \(0\le k < N-2\) . For \(\mu \le (\frac{N-k-2}{2})^2\) , put \(L_\mu := \Delta + \mu d_{\Sigma }^{-2}\) where \(d_{\Sigma }(x) = \textrm{dist}(x,\Sigma )\) . We study boundary value problems for equation \(-L_\mu u = g(u,|\nabla u|)\) in \(\Omega \setminus \Sigma \) , subject to the boundary condition \(u=\nu \) on \(\partial \Omega \cup \Sigma \) , where \(g: {\mathbb {R}} \times {\mathbb {R}}_+ \rightarrow {\mathbb {R}}_+\) is a continuous and nondecreasing function with \(g(0,0)=0\) , \(\nu \) is a given nonnegative measure on \(\partial \Omega \cup \Sigma \) . When g satisfies a so-called subcritical integral condition, we establish an existence result for the problem under a smallness assumption on \(\nu \) . If \(g(u,|\nabla u|) = |u|^p|\nabla u|^q\) , there are ranges of p, q, called subcritical ranges, for which the subcritical integral condition is satisfied, hence the problem admits a solution. Beyond these ranges, where the subcritical integral condition may be violated, we establish various criteria on \(\nu \) for the existence of a solution to the problem expressed in terms of appropriate Bessel capacities.