Unified lattice sums accommodating multiple sublattices for solutions of the Helmholtz equation in two and three dimensions

The evaluation of the interaction between objects arranged on a lattice requires the computation of lattice sums. A scenario frequently encountered involves systems governed by the Helmholtz equation in the context of electromagnetic scattering in an array of particles forming a metamaterial, a metasurface, or a photonic crystal. While the convergence of direct lattice sums for such translation coefficients is notoriously slow, the application of Ewald's method converts the direct sums into exponentially convergent series. We present a derivation of such series for the two-dimensional (2D) and three-dimensional (3D) solutions of the Helmholtz equation, namely, spherical and cylindrical solutions. When compared to prior research, our expressions are especially aimed at computing the lattice sums for several interacting sublattices in one-dimensional lattices (chains), 2D lattices (gratings), and 3D lattices. The presented approach establishes a unified framework to derive efficient expressions for lattice sums. Besides reproducing previously known formulas in a more accessible manner the approach results in expressions for interacting sublattices of a chain in three dimensions. Furthermore, the derived formulas are not limited to the dipolar case but are applicable to arbitrary multipolar orders. We verify our results by comparison with the direct computation of the lattice sums.

Figure 1

Illustrations for contemporary photonic materials featuring multiple sublattices. From panels (a) to (d) these are zigzag chains that mimic the Su-Schrieffer-Heeger model, metasurfaces where bound states in continuum can be observed, plasmonically induced transparency in dolmen structures, and moiré lattices.

Figure 2

Layout of the geometry in the lattice summations. In all cases, we show two different sublattices with blue dots and green stars. The shift between these sublattices is described by $\mathit{r}$. The dotted lines indicate the decomposition of that shift into ${\mathit{r}}_{\perp }$ and ${\mathit{r}}_{\parallel }$. In panels (a), (b), and (d), the different sublattices are not required to be in the same plane or line. The vector $\mathit{R}$ is a lattice vector. The first row with panels (a) to (c) shows the case of spatial dimension $d=3$ and $d^{\prime }\in \{1,2,3\}$, respectively. In all three cases we use the coordinate system indicated on the left. The second row with panels (d) and (e) shows the cases for $d=2$, again, with coordinates as shown on the left. The orange color shows the spatial domain of the lattice summation. The parallel component of the wave vector ${\mathit{k}}_{\parallel }$ must lie in this domain.

Figure 3

Comparison of the direct evaluation of the series and the value for the exponentially convergent expressions. Each panel shows the real part (solid line) and imaginary part (dashed line) of the value for the direct summation (blue), the exponentially convergent ones (orange), and the relative deviation of the direct summation (green). The $x$ axis shows the number of included layers. These layers have a square and cubic shape for 2D and 3D lattices, respectively. Panels (a), (b), and (c) show the values for the 1D lattice (chain), 2D lattice (grating), and 3D lattice examples for the spherical solution of the Helmholtz equation, respectively. Panels (d) and (e) show the values for the 1D lattice and 2D lattice for the cylindrical solution, respectively. Panel (f) shows an example of an application of the lattice sums for a chain of spheres on two sublattices. The lattice sums are used to calculate the coupling of the spheres within the T-matrix framework and also to compute the electric field the intensity of which is shown here

Figure 4

Comparison of the direct summation and the value of exponentially fast converging series for large shifts ${\mathit{r}}_{\perp }$ perpendicular to the lattice. Panel (a) and (b) show the results for 1D and 2D lattice summing spherical solutions. The values of $l=2$ and $m=1$ are chosen, which in case of ${\mathit{r}}_{\perp }=0$ must vanish. Panel (c) shows the case $d=2$ and $d^{\prime }=1$.

Figure 5

The same parameters as in Fig. 3 are used but here we use a convolution to average over the oscillations to obtain a faster convergence of the direct summation. Still, the direct approach needs significantly longer.