What are good non-paraxial gaussian-beam-like solutions of the Helmholtz equation?

Question

I am playing around with some optics manipulations and I am looking for beams of light which are roughly gaussian in nature but which go beyond the paraxial regime and which include non-paraxial vector-optics effects like longitudinal polarizations and the like.

That is, I am looking for monochromatic solutions of the Maxwell equations which look roughly like a tightly-focused gaussian beam and which

have finite energy in each transversal slab,
have analytical forms that are manageable enough to work with,
are full solutions to the Maxwell equations (i.e. not paraxial, or some subleading-order approximation to it),
fully include vector effects, all the way through to the forwards-circular polarization states mentioned e.g. here, and
ideally, come as a full basis that's separated along orbital angular momentum, as the Laguerre-Gauss solutions do.

In the past I've used leading-order corrections to the paraxial approximation (like the ones used here, which were drawn from here), and those work fine when the non-paraxiality is small and the forwards ellipticity can be considered perturbative, but I would like to take things to situations where the non-paraxial forwards ellipticity is large enough to induce circular polarizations while still maintaining a manageable analytical treatment that's also fully rigorous.

Is this doable? Or must the passage to a full non-paraxial solution require the use of infinite-energy Bessel-like beams if one wants to maintain analytical solvability?

Pointers to references with detailed work-throughs of such solutions would be highly appreciated - I want to use these as tools to build stuff on top of, and I'd like a solid base to work on.

Emilio Pisanty · Accepted Answer

The cleanest way to do this is using what is known as Complex Focus fields, which are built using a small set of simple (though non-obvious) key ideas:

The basic ingredient is the multipolar solution of the Helmholtz equation,
$$Lambda_{l,m}(mathbf r) = 4pi i^l : j_l(kr)Y_{lm}(theta,phi),$$
where $j_l(kr)$ is a spherical Bessel function and $Y_{lm}(theta,phi)$ is a spherical harmonic, and which satisfies $(nabla^2+k^2)Lambda=0$.
In terms of the Cartesian coordinates of $mathbf r=(x,y,z)$, this looks somewhat intimidating, as $theta$ and $phi$ are discontinuous functions and $j_l(kr)$ has a square root, but we also know that $j_l(kr)sim r^l$ times a convergent series in $r^2$, and if this factor is incorporated into the spherical harmonic we get the solid spherical harmonic $S_{lm}(mathbf r) = r^l : Y_{lm}(theta,phi)$, which is a homogeneous polynomial of degree $l$ in $x,y,z$.
... all of which is a long-winded way of saying that $Lambda_{l,m}(mathbf r)$ is not just a continuous function of $x,y,z$, but it is actually an entire function.

The second key ingredient is using the fact that $Lambda_{l,m}(mathbf r)$ is entire, and displacing the $z$ coordinate by an imaginary offset, to $Lambda_{l,m}(mathbf r-izetahat{mathbf{e}}_z)$, which does not affect the Helmholtz equation $(nabla^2+k^2)Lambda=0$.
As $zeta$ increases, the solutions morph away from spherical-wave character, and their angular spectrum becomes more and more concentrated on forwards-propagating waves. In the large $kzeta$ limit, the solutions first become tightly-focused nonparaxial beams, and then they de-focus to approach propagating beams in the paraxial limit. More precisely,

the monopolar solution $Lambda_{0,0}(mathbf r-izetahat{mathbf{e}}_z)$ approaches a Gaussian beam, and
the 'extremal' cases $Lambda_{l,pm l}(mathbf r-izetahat{mathbf{e}}_z)$ with nonzero angular momentum approach the (base) Laguerre-Gauss beams without radial nodes.
(The middle cases, $Lambda_{l,m}(mathbf r-izetahat{mathbf{e}}_z)$ with $|m|<l$, become complex standing-wave mixtures of waves in various directions.)

Finally, to produce suitable vector solutions (as opposed to the scalar $Lambda_{l,m}$) which satisfy the transversality condition $nablacdotmathbf E=0$ in addition to the Helmholtz equation, one can use suitable differential operators (with good examples being $mathbf V_{mathbf p} f(mathbf r) = tfrac{1}{ik} nabla times (mathbf p f(mathbf r))$ and $mathbf V_{mathbf p} f(mathbf r) = -tfrac{1}{k^2} nabla times left(nabla times (mathbf p f(mathbf r))right)$) which involve single derivatives and as such do not affect the Helmholtz equation.

A good reference that delves into these fields in depth is

Scalar and electromagnetic nonparaxial bases composed as superpositions of simple vortex fields with complex foci. R. Gutiérrez-Cuevas and M.A. Alonso. Opt. Express 25, 14856 (2017).

and my own usage will come out shortly as

Optical polarization skyrmionic fields in free space. R. Gutiérrez-Cuevas and E. Pisanty. Under review (2020).

My implementation on Mathematica is available as the ComplexFocus package on GitHub, over at github.com/ComplexFocus/ComplexFocus.

akhmeteli · Answer

I would think Section II of my old work,

Efficient Heating of Thin Cylindrical Targets by Broad Electromagnetic Beams I. Andrey Akhmeteli. arXiv:physics/0405091 (2004).

contains exactly what you need: a relatively simple exact solution of the free Maxwell equations that is asymptotically approximated by a Gaussian beam in the limit $delta/lambdarightarrowinfty$, where $delta$ is the size of the waist of the Gaussian beam and $lambda$ is the wavelength.
The solution is expressed in terms of Hertz potentials (see formulas 1-6, 21-22). The solution is written for a circularly polarized Gaussian beam, but it is not difficult to modify the solution to get a linearly polarized beam.
This solution has the form of a single 1D integral, which must be integrated numerically (as it seems that Mathematica cannot integrate it symbolically). The paper presents the Hertz potentials as an integral; if you want the fields directly, you should take the derivatives in the double curl symbolically before the numerical integration, and you should be careful to use formulas for derivatives of the Bessel functions which do not cause loss of precision (so that you do not compute a difference of two functions which do not differ much in the vicinity of zero).
The numerical integral is oscillatory, but it is up to you to decide whether it is "highly oscillatory". You can compute it using the fast Fourier transform or the fast Hankel transform (Fourier-Bessel transform), depending on whether you need values along a line parallel to the axis of the beam, or along a radius.

What are good non-paraxial gaussian-beam-like solutions of the Helmholtz equation?

2 Answers

Add your own answers!

Ask a Question