  We present a solver for the 2D high-frequency Helmholtz equation in
heterogeneous, constant density, acoustic media, with online parallel
complexity that scales empirically as $\mathcal{O}(\frac{N}{P})$, where $N$ is
the number of volume unknowns, and $P$ is the number of processors, as long as
$P = \mathcal{O}(N^{1/5})$. This sublinear scaling is achieved by domain
decomposition, not distributed linear algebra, and improves on the $P
=\mathcal{O}(N^{1/8})$ scaling reported earlier in [L. Zepeda-N\'u\~nez and L.
Demanet, J. Comput. Phys., 308 (2016), pp. 347-388 ]. The solver relies on a
two-level nested domain decomposition: a layered partition on the outer level,
and a further decomposition of each layer in cells at the inner level. The
Helmholtz equation is reduced to a surface integral equation (SIE) posed at the
interfaces between layers, efficiently solved via a nested version of the
polarized traces preconditioner [L. Zepeda-N\'u\~nez and L. Demanet, J. Comput.
Phys., 308 (2016), pp. 347-388.]. The favorable complexity is achieved via an
efficient application of the integral operators involved in the SIE.
