  We present a novel computational methodology for solving the scalar nonlinear
Helmholtz equation (NLH) that governs the propagation of laser light in Kerr
dielectrics. The methodology addresses two well-known challenges in nonlinear
optics: Singular behavior of solutions when the scattering in the medium is
assumed predominantly forward (paraxial regime), and the presence of
discontinuities in the % linear and nonlinear optical properties of the medium.
Specifically, we consider a slab of nonlinear material which may be grated in
the direction of propagation and which is immersed in a linear medium as a
whole. The key components of the methodology are a semi-compact high-order
finite-difference scheme that maintains accuracy across the discontinuities and
enables sub-wavelength resolution on large domains at a tolerable cost, a
nonlocal two-way artificial boundary condition (ABC) that simultaneously
facilitates the reflectionless propagation of the outgoing waves and forward
propagation of the given incoming waves, and a nonlinear solver based on
Newton's method.
  The proposed methodology combines and substantially extends the capabilities
of our previous techniques built for 1Dand for multi-D. It facilitates a direct
numerical study of nonparaxial propagation and goes well beyond the approaches
in the literature based on the "augmented" paraxial models. In particular, it
provides the first ever evidence that the singularity of the solution indeed
disappears in the scalar NLH model that includes the nonparaxial effects. It
also enables simulation of the wavelength-width spatial solitons, as well as of
the counter-propagating solitons.
