  We study the nonlinear Helmholtz equation $(\Delta - \lambda^2)u = \pm
|u|^{p-1}u$ on $\mathbb{R}^n$, $\lambda > 0$, $p \in \mathbb{N}$ odd, and more
generally $(\Delta_g + V - \lambda^2)u = N[u]$, where $\Delta_g$ is the
(positive) Laplace-Beltrami operator on an asymptotically Euclidean or conic
manifold, $V$ is a short range potential, and $N[u]$ is a more general
polynomial nonlinearity. Under the conditions $(p-1)(n-1) > 4$ and $k >
(n-1)/2$, for every $f \in H^k(S^{n-1}_\omega)$ of sufficiently small norm, we
show there is a nonlinear Helmholtz eigenfunction taking the form
\begin{equation*} u(r, \omega) = r^{-(n-1)/2} \Big( e^{-i\lambda r} f(\omega) +
e^{+i\lambda r} b(\omega) + O(r^{-\epsilon}) \Big), \qquad \text{as } r \to
\infty, \end{equation*} for some $b \in H^k(S_\omega^{n-1})$ and $\epsilon >
0$. That is, the scattering matrix $f \mapsto b$ preserves Sobolev regularity,
which is an improvement over the authors' previous work with Zhang, that proved
a similar result with a loss of four derivatives.
