1. Introduction

Interactions between structured optical fields (SOFs) and scatterers produce many interesting effects including optical manipulations [1], selective excitations of multipolar resonances [2–6], angular momentum (AM) dichroisms [7–10], and breaking of the conventional optical theorem [11]. Examples of such SOFs include localized beams with non-zero longitudinal field components [12], vortex beams with intrinsic orbital angular momentum (OAM) and helical phasefronts [13], and cylindrical vector (CV) beams with spatially inhomogeneous polarization states [14].

Recently, there have been advancements on the generation of SOFs using spatial light modulators and metasurfaces [15], and on the synthesis of particles with peculiar morphologies [16] and highly anisotropic responses. These advancements on meta-photonics are expected to offer unprecedented light matter interactions that could not be realized using planewaves or nonresonant Rayleigh scatterers. For instance, AM dichroisms involving meta-atoms and vortex beams have been recently exploited for excitation of dark modes [9], enantioselective sensing of a chiral dipolar helix using vortex beams with OAM [8], and strong OAM dichroism of meta-atoms [10]. Despite these possibilities, it is difficult to study light matter interactions involving both SOFs and arbitrarily shaped particles, and especially research involving nonspherical particles are still at its infancy.

Accurate and efficient description of both shaped beams and arbitrary scatterers is strongly required to study this problem, and multipole expansion (ME) can be a versatile framework for this task. Under this framework, a scatterer is described as a $T$-matrix [17], which can be calculated for arbitrarily shaped particles and very efficiently for particles with rotational symmetry [18]. Because only a few multipole order is enough to accurately describe the scattering phenomena involved with subwavelength particles [19,20], the subwavelength particles are efficiently approximated as low-order multipoles under ME. An arbitrarily shaped beam is incorporated into the framework as the beam-shape coefficients (BSCs) [21,22]. Interaction between arbitrarily shaped beams and spherical particles has been studied under the generalized Lorenz-Mie theory [23].

In this work, we first summarize some SOFs, where explicit expressions of their fields and BSCs are provided and their properties are briefly discussed. The discussed SOFs include Laguerre-Gaussian (LG) beams, tightly-focused LG beams, and Bessel beams, and CV beams. Then, selective excitations of multipolar resonances of a sphere using vortex beams with AM is discussed using ME. Also, AM dichroisms involving a chiral sphere and a twisted-nanorods are discussed using ME, where multipolar resonances of desired order and parity, as well as orientation, can be selectively excited using engineered SOFs. We hope this study can be helpful for efficiently studying the interaction between SOFs and structured meta-atoms.

2. Structured optical fields and their beam-shape coefficients

2.1 Laguerre-Gaussian beams

First, we recall the description of propagating SOFs. One of the simplest and widely adopted beam models would be the paraxial LG beam, whose electric field is given as [13]

The paraxial expression (Eq. (1)) provides a transparent visualization of shaped beam structures. However, the paraxial beams are essentially approximate solutions of the Maxwell equations and cannot be multipole decomposed; that is, they cannot be expressed as the sum of basis functions that satisfy the Maxwell equations (Eq. (13)). A Maxwellian LG beam can be obtained using the far field matching [21,24] as

This Maxwellian LG beam can be multipole decomposed by inserting Eq. (2) into Eq. (15). Using the integral representation of the Bessel function $\int _0^{2\pi }e^{in\beta }e^{ix\cos {\beta }}\mathrm {d}\beta = 2\pi i^n J_n(x)$, BSC around an arbitrary origin $\mathbf {r}_p$ can be efficiently calculated [22]. Throughout this work, $(\rho _p,\phi _p,z_p) = \mathbf {r}_p$ in cylindrical coordinate assuming the beam center is at the origin. The resulting BSC is

It is worthwhile to note that Eq. (2c) at $p = 0$ and $l = 0$ exactly recovers the nonparaxial Gaussian beam obtained using the vector angular spectrum representation (ASR) ignoring the evanescent field contribution given as [25]

The fields obtained using the vector ASR and their BSCs can be found in Ref. [23], but this method may be useful for a limited beams with known angular spectra, such as simple Gaussian beams [26], Bessel-Gaussian beam [27], LG beam, and Airy beam [28]. Collimated beams with only tangential components of the fields are often discussed using the ASR, but they do not satisfy the Maxwell equations and cannot be exactly multipole decomposed. The obtained BSCs using this method are different from Eq. (5) with $\cos {\alpha }~\pi _{nm}(\alpha )$ and $\cos {\alpha }~\tau _{nm}(\alpha )$ changed into $\pi _{nm}(\alpha )$ and $\tau _{nm}(\alpha )$, respectively [26,27].

Unlike planewaves, a light beam has a finite beam width and generally diffracts or diverges away from its focal point. Near the laterally and axially localized focal point, objects can be trapped using the gradient force [1]. For simple Gaussian beam with $p = 0$ and $l = 0$ (Fig. 1(a–b)), beam divergence angle is given as $\alpha _0 = \arctan [2/(kw_0)]$, which relates to the numerical aperture (NA) as $\mathrm {NA} = n\sin {\alpha _0}$. A simple Gaussian beam with $w_0 = 10\lambda _0$ has $\alpha _0 = 1.8^\circ$ with a collimated beam profile (Fig. 1(a)), and a more localized Gaussian beam with $w_0 = \lambda _0$ has $\alpha _0 = 17.7^\circ$ with a diverging beam profile (Fig. 1(b)).

 

Fig. 1. Laguerre-Gaussian beams. Electric field intensity $|\mathbf {E}|^2$ and polarization states $(E_x,E_y)$ (arrows) on $xy$-plane at $z=0$ (top) and $xz$-plane at $y=0$ (bottom) calculated using Eq. (2c). Blue and green arrows have $\pi /2$ phase difference. Insets are phase of $E_x$. Beam parameters are: $\lambda _0$ = 1 $\mu$m. $(p_x,p_y) = (1,0)$.

Download Full Size | PDF

In the paraxial limit, a LG beam with $l\neq 0$ is a vortex beam with well-defined intrinsic OAM [13], which can be transferred to objects producing an orbital motion with respect to the beam axis [1,29]. The characteristic rotating phasefronts have the phase singularity at the beam axis with the rotation direction depending on the sign of $l$ (insets of Fig. 1(c–e)). Unlike spin angular momentum (SAM) of photons $s\hbar$ with two possible values $s=\pm 1$, OAM $l\hbar$ has a large degree-of-freedom with $l = \cdots ,-2,-1,0,+1,+2,\ldots$, which can be applied to communications, encryption [30], and holography [31]. Vortex beams have been studied for OAM dependent light matter interactions for characterization of higher-order chiral molecular transitions [7] and strong AM dichroisms of chiral meta-atoms [8–10]. In addition, vortex beams with annular beam intensity profiles (Fig. 1(c–e)) can be used for particle trapping [1] and stimulated emission depletion microscopy [32]. Another degree-of-freedom radial index $p$ creates radial fringes with $p+1$ concentric rings and a reduced focal beam size (Fig. 1(f)).

However, it should be noted that nonparaxial beams discussed in this work do not have well-defined orbital and spin AMs [33], because the operators for SAM and OAM break the transversality of Maxwell fields [34]. In this work, $s=\pm 1$ denote $(p_x,p_y)=(1,\pm i)/\sqrt {2}$ that correspond to circularly-polarized states in the paraxial limit, and $s=0$ denotes $(p_x,p_y)=(1,0)$.

2.2 Tightly-focused Laguerre-Gaussian beams

The Debye-Wolf vector diffraction theory has been popularly used to obtain Maxwellian beams tightly-focused by aplanatic lens, whose focal fields are given as [35,36]

The complex polarization vector $\hat {\mathbf {Q}}$ is the polarization of the field after passing through the focusing lens and is given as

 

Fig. 2. Tightly-focused Laguerre-Gaussian beams. Electric field intensity $|\mathbf {E}|^2$ and polarization states $(E_x,E_y)$ (arrows) on $xy$-plane at $z=0$ (top) and $xz$-plane at $y=0$ (bottom) calculated using Eq. (9a). Blue and green arrows have $\pi /2$ phase difference. Inset is phase of $E_x$. Beam parameters are: $\lambda _0$ = 1 $\mu$m, $p = 0$, $f$ = 1 mm, $w$ = 1 mm.

Download Full Size | PDF

Compared to LG beams (Eq. (2)), tightly-focused LG beams (Eq. (7)) are circularly-symmetric in terms of the Poynting vector and energy density distributions [37]. Electric and magnetic fields of LG beam have quite different functional form (Eq. (2)), compared to those of tightly-focused LG beam (Eq. (7)). In addition, Eq. (2c) results $E_x = 0$ at $p_y = 0$ and $E_y = 0$ at $p_x = 0$, whereas this is not true for Eq. (7a), as well as magnetic fields (Eq. (2b) and (7b)).

2.3 Circularly-symmetric Bessel beams

Bessel beams have been a versatile theoretical beam model for radially localized beams (Fig. 3(a–b)), vortex beams, and CV beams (Fig. 3(c–f)), because Bessel beams are Maxwellian with analytical expressions (Eq. (9) and (11)). Bessel beams are also useful for scattering calculations, because their BSCs are also analytically expressed (Eqs. (10) and (12)).

 

Fig. 3. Bessel beams and CV beams. Electric field intensity $|\mathbf {E}|^2$ and polarization states $(E_x,E_y)$ (arrows) on $xy$-plane at $z=0$ (top) and $xz$-plane at $y=0$ (bottom) calculated using Eq. (9a) and (11a). Beam parameters are: $\lambda _0$ = 1 $\mu$m, $l = 0$; (c–f) $\alpha _0 = 30^\circ$.

Download Full Size | PDF

Using the angular spectrum $\hat {\mathbf {E}}(\alpha ,\beta )=\tfrac {\delta (\alpha -\alpha _0)}{\sin {\alpha _0}}e^{il\beta }\hat {\mathbf {Q}}(\alpha ,\beta )$ [38], the fields of circularly-symmetric Bessel beams are expressed as

Maxwellian Bessel beams with simpler forms have been obtained based on the vector potential [9,10], resulting functional forms similar to LG beams discussed before (Eq. (2)). Circularly symmetric Bessel beams have been obtained by using the vector ASR or by averaging dual fields obtained using the vector potential [37].

Bessel beams have been noticed as diffraction-less beams [38,39]; this feature can be confirmed in Eq. (9) and (11) where the expressions have argument $z$ separable as $e^{ik\cos \alpha _0z}$, which can also be confirmed in Fig. 3 where the field profiles do not change as they propagate in $z$ direction. Bessel beams self-heal even after encountering obstacles; this feature allows simultaneous optical trapping of array of particles over a large distances [39].

2.4 Cylindrical vector beams

A CV beam, or cylindrically polarized beam, has a spatially inhomogeneous polarization profile given as $(p_x,p_y) = p_\rho (\cos \beta ,\sin \beta )+p_\phi (-\sin \beta ,\cos \beta )$, where $p_\rho$ and $p_\phi$ determine the polarization state of the CV beam. For special cases, $p_\phi =0$ corresponds to a radially polarized beam with $(p_x,p_y) = (\cos \beta ,\sin \beta )$, and $p_\rho =0$ corresponds to an azimuthally polarized beam with $(p_x,p_y) = (-\sin \beta ,\cos \beta )$ [41]. Noting that $(\cos \beta ,\sin \beta )=e^{+i\beta }(1,-i)/2+e^{-i\beta }(1,+i)/2$ and $(-\sin \beta ,\cos \beta )=ie^{+i\beta }(1,-i)/2-ie^{-i\beta }(1,+i)/2$, a CV beam can be interpreted as a superposition of two vortex beams with $l=+1, s=-1$ and $l=-1, s=+1$. It is worthwhile to note that vector beams with $l=+1$ and $s=-1$ consists of out-of-phase radial and azimuthal polarizations (Fig. 2(f)).

For simplicity, we express CV beams using the circularly-symmetric Bessel beams as

Because a CV beam consists of vortex beams, it exhibits doughnut-shaped beam profile even with $l = 0$ (Fig. 3(c–f)). Throughout this work, we will only consider this $l = 0$ case for CV beams. An important feature of CV beams is their well-defined parity. It is straghtforward to see from Eq. (11) that a radially polarized beam is TM mode with $H_z = 0$, and an azimuthally polarized beam is TE mode with $E_z = 0$ [42]. Because of this longitudinal field, a radially polarized beam (Fig. 3(c)) has a nonzero field intensity at the beam axis, whereas an azimuthally polarized beam (Fig. 3(d)) has a zero field intensity at the beam axis. Tightly-focused radially polarized beams exhibit sharp beam spots below the Abbe’s diffraction limit for longitudinal field components [43] for super-resolution imaging. Simple linear combinations of radially and azimuthally polarized states can yield spirally polarized states with spirally rotating polarizations (Fig. 3(e–f)).

3. Interactions between structured optical fields and meta-atoms

Multipolar resonances are at the heart of the scattering effects from meta-atoms. Recently investigated structured meta-atoms have shown their potentials in controlling multipolar resonances at planewave incidences. Notably, coupled-plasmonic and high-refractive-index subwavelength particles with strong higher-order multipolar resonances have been engineered to exhibit interesting optical effects including directional scattering, optical magnetism, strong artificial chirality, and nonradiative anapoles [20,44]. Another method to control the multipolar resonances is to engineer the incident light beams, which we discuss in this section [2–6]. Incorporating both approaches, generalized cases involving complex SOFs and nonspherical meta-atoms would potentially exhibit unprecedented optical effects.

3.1 Selective excitations of multipolar resonances

First, we discuss selective excitations of multipolar resonances of a sphere using vortex beams with AM. Scattering of a sphere by tightly-focused beams with AM exhibits excitation of fewer multipolar resonances [4] with spectrally sharper scattering peaks [5]. In typical Mie scatterings, a high-refractive-index dielectric sphere exhibits distinct multipolar resonances, which can be observed at the planewave incidence (Fig. 4(b)). Unlike this planewave incidence that excite all multipolar resonances, vortex beams with $\pm |l|$ and $s=\pm 1$ excites only the multipolar modes of a sphere higher than $|l|+1$ (Fig. 4(c)); for instance, $l=+2$ and $s=+1$ case does not excite dipolar and quadrupolar resonances of the dielectric sphere (red line of Fig. 4(c)). By inhibiting low-order multipolar resonances using vortex beam with AM, sharper and often weaker higher-order multipolar resonances can be selectively excited [3].

In addition, an electric or magnetic mode can be selectively excited using a light beam with radial or azimuthal polarization, respectively (Fig. 4(d)) [2,3]. For scattering by CV beams, conventional optical theorem may need to be revised [11]. Because radially polarized beam has only longitudinal electric field component at the beam axis (Eq. (11a)), a small plasmonic sphere excited by this beam would have only an electric dipole mode in the longitudinal direction, which will radiate into the transverse directions.

 

Fig. 4. Selective excitations of multipolar resonances of a dielectric sphere using vortex Bessel beams. (a) Schematics of a dielectric sphere illuminated by vortex Bessel beam with NA = 0.8. (b) Multipole decomposed and total scattering of the sphere by planewave incidence. (c, d) Scattering of the sphere illuminated by (c) vortex beams and (d) CV beams.

Download Full Size | PDF

This selective excitation of multipolar resonances has been used to control the polarization of third harmonic generation [45]. Another application is the excitation of ideal nonradiating anapole, where counter-propagating azimuthally polarized beams could only excite a single multipolar resonance with one parity [6]. In addition, SOFs with AM can excite of dark modes, which are otherwise not accessible by planewaves [9].

3.2 Angular momentum dichroisms

Complementary chiral objects, or enantiomers, have been typically distinguished using circular dichroism, or the SAM dependent absorption of a system due to opposite circularly-polarized planewaves. This SAM dichroism can be extended using vortex beams that carry OAM. To study AM dichroisms, we utilize four different dichroisms defined in Ref. [10] as: SAM dichroism $\Delta \sigma ^\mathrm {SAM}_{l} = \sigma _{l,+1}-\sigma _{l,-1}$, OAM dichroism $\Delta \sigma ^\mathrm {OAM}_{l,s} = \sigma _{+|l|,s}-\sigma _{-|l|,s}$, parallel class dichroism $\Delta \sigma ^\mathrm {PC}_{l} = \sigma _{+|l|,+1}-\sigma _{-|l|,-1}$, and anti-parallel class dichroism $\Delta \sigma ^\mathrm {APC}_{l} = \sigma _{+|l|,-1}-\sigma _{-|l|,+1}$, where $\sigma _{l,s}$ is extinction from a vortex beam with $l$ and $s$. Note that spin and orbital AM are not well-defined in the nonparaxial beams as mentioned before; nevertheless, we use the terms OAM and SAM dichroisms following the convention [10]. Because vortex beams are mirror symmetric if signs of both $l$ and $s$ are reversed, $\Delta \sigma ^\mathrm {SAM}_l$ and $\Delta \sigma ^\mathrm {OAM}_s$ cannot be used for chiral sensing if $l\neq =0$ and $s\neq 0$, respectively. For instance, vortex beams with $l\neq 0$ and $s=\pm 1$ are not mirror symmetric to each other, so they can produce nonzero $\Delta \sigma ^\mathrm {SAM}_l$ for a non-chiral object [46].To study AM dichroisms using vortex beams, we utilize a lossless weakly chiral sphere with $n=3.5$ and chirality parameter $\kappa =0.01$. The AM dichroism spectra (Fig. 5(b)) consist of peaks at the multipolar resonances of the reference dielectric sphere (Fig. 4(b)). $\Delta \sigma ^\mathrm {SAM}_0$ is relevant to the conventional circular dichroism at the paraxial limit, and is strong near 400, 550, and 850 THz (black line of Fig. 5(b)) where dipolar resonances exist (Fig. 4(b)). $\Delta \sigma ^\mathrm {PC}_1$ inhibits the differential scattering originating from dipolar resonances (red line of Fig. 5(b)), because $\sigma _{\pm |l|,\pm 1}$ do not excite dipolar resonances. Therefore, $\Delta \sigma ^\mathrm {PC}_1$ can be used to observe chiroptical responses related to quadrupolar resonances by suppressing the dipolar responses, which are often broader and stronger. $\Delta \sigma ^\mathrm {OAM}_0$ has been used to observe chiral light matter interactions involving quadrupolar responses [7], but it can also be used to observe chiroptical responses from dipolar chiral scatterer [8]. Also, $\Delta \sigma ^\mathrm {APC}_1$ can access both dipolar and higher-order multipolar resonances.

 

Fig. 5. (a) Schematics of a weakly chiral dielectric sphere illuminated by Bessel beam. (b) AM dichroisms $\Delta \sigma ^\mathrm {SAM}_0$, $\Delta \sigma ^\mathrm {OAM}_0$, $\Delta \sigma ^\mathrm {PC}_1$, and $\Delta \sigma ^\mathrm {APC}_1$ of the chiral sphere.

Download Full Size | PDF

This problem of chiral interactions between dipolar chiral objects and vortex beams with $|l|>0$ can be examined using the local optical chirality density of the fields [8]. Microscopically, enantioselective extinction of dipolar chiral objects depends on the local optical chirality density of the fields [47]. This also applies to general SOFs. Here, we separated the optical chirality density into transverse and longitudinal origins as: $S_{xy}=\eta \mathrm {Im}(E_xH_x^*+E_yH_y^*)/E_0^2$ and $S_{z}=\eta \mathrm {Im}(E_zH_z^*)/E_0^2$ (Fig. 6). The optical chirality density was normalized by that of planewaves for simplicity. A planewave propagating in $z$ direction would only have $S_{xy}$, but tightly-focused beams may have strong $S_{z}$.

 

Fig. 6. Normalized optical chirality densities in transverse and longitudinal parts: $S_{xy}$ and $S_{z}$ of Bessel beams with NA = 0.9 calculated using Eq. (9).

Download Full Size | PDF

A Bessel beam with $l=0$ and $s=+1$ (Fig. 6(a)) has strong $S_{xy}$ at the focal point, which relates to the conventional circular dichroism in the paraxial limit. Interestingly, a vortex Bessel beam with $l=+1$ and $s=0$ (Fig. 6(b)) has $S_{z}$ at the focal point. This attributes to OAM-to-SAM conversion due to tight-focusing of non-paraxial beam [29,33], which results non-zero $\Delta \sigma ^\mathrm {OAM}_0$ at dipolar resonances of the chiral sphere (Fig. 5(b)) and dipolar helix [8].

The tight-focusing is responsible for the spin-orbit interaction where part of SAM is converted to OAM [29,33], and the reverse is also possible [8]. The longitudinal field components demonstrate this spin-orbit interactions as the doughnut-like field pattern from simple circularly-polarized beams (Fig. 6(a)) and spin-polarized focal fields from linearly-polarized vortex beams (Fig. 6(b)). The magnitude of the longitudinal fields can be controlled by the degree of focusing [12]. In addition, the higher-order harmonic generation generated from longitudinal field components can have different polarization states compared to transverse field components [48].

Notably, vortex beams with $l = +1$ and $s = +1$ have zero optical chirality density at the focal point (Fig. 6(c)), resulting no AM dichroism originating from the dipolar resonances (Fig. 5(b)). On the other hand, vortex beams with $l = +1$ and $s = -1$ have strong optical chirality density in the longitudinal direction (Fig. 6(d)), and $\Delta \sigma ^\mathrm {APC}_1$ appear in the dipolar resonances (Fig. 5(b)).

3.3 Angular momentum dichroisms of chiral meta-atoms with anisotropic responses

Structured meta-atoms often have highly anisotropic responses. A subwavelength plasmonic twisted-nanorods [10,44] have dipolar resonances anisotropically in $x$ and $y$ directions, which correspond to the broader resonance at 560 nm and the sharper resonance at 570 nm, respectively (Fig. 7(a)).

 

Fig. 7. Extinction $\sigma$ and differential extinction $\Delta \sigma$ of a twisted-nanorods illuminated by Bessel beams in different directions denoted by different colors. (a) Schematics. (b–d) $\sigma$ and $\Delta \sigma$ by (b) circularly-polaried Bessel beams with $l=0$ and $s=\pm 1$, (c) vortex Bessel beams with $l=\pm 1$ and $s=\pm 1$, and (d) vortex Bessel beams with $l=\pm 1$ and $s=\mp 1$. $\sigma$ and $\Delta \sigma$ in (c) are enlarged by 10 and 50 times for better visualization. NA = 0.9.

Download Full Size | PDF

In general, a planewave excites transverse modes (e.g., $x$ and $y$), which may not be degenerate for nonspherical meta-atoms. Using SOFs, we can selectively excite only longitudinal mode where the transverse modes are inhibited. We calculated extinction at three different relative orientations between the beam axis and the meta-atom, and the orientations are denoted by different colors (Fig. 7).

Under the illumination of Bessel beams with $l = 0$ and $s = \pm 1$, transverse modes are excited (Fig. 7(b)). Beams propagating in $z$-direction (black lines) excite both $x$ and $y$ resonances, whereas beams propagating in $x$- or $y$-direction (red or blue lines) excite either $y$ or $x$ resonance.

Vortex beams with $l=\pm 1$ and $s=\pm 1$ cannot directly excite the dipolar resonances (Fig. 6(c)), but weak resonances at the dipolar resonances are observed (Fig. 7(c)). These responses are excited through quadrupolar fields [49]; this transition is weak, but allowed for chiral particles [44]. Finally, vortex beams with $l=\pm 1$ and $s=\mp 1$ directly excite the dipolar chiroptical responses in the longitudinal direction (Fig. 6(d)). Because the responses parallel to the beam axis are excited, so the beams propagating in $z$-direction do not see the resonances (Fig. 7(d)). Therefore, we can selectively excite transverse or longitudinal responses of meta-atoms using SOFs.

4. Conclusions

In this work, we first presented a summary of some SOFs and their BSCs; the discussed beams include LG, tightly-focused LG, and Bessel beams, and CV beams. The discussed beams are Maxwellian and can be multipole decomposed [50]. The field expressions presented in this work can be useful by themselves, because fields can be directly inserted into some fully numerical methods as the background (input) fields. For some other methods, such as FDTD, fields may be inserted as a discrete sum of planewaves whose weight is given by the angular spectrum.

In addition, selective excitation of multipolar resonances using SOFs with AM was discussed. Multipolar resonances of the particular order and parity could be excited using the engineered SOFs, and the dipoles in the longitudinal direction can be selectively excited using SOFs. Also, AM dichroisms arising from the interaction between chiral meta-atoms and SOFs were discussed to illustrate ME as the potential theoretical tool for studying meta-photonics involving SOFs.

In addition, this work can also be helpful for metasurfaces. Typically, ME has been used to describe a single particle or small clusters [20], but ME can describe electromagnetically interacting particles dispersed over a large physical space in a rigorous and efficient manner beyond the single particle or cluster level. This feature makes ME potentially useful for studying metasurfaces and particulate media. ME can rigorously describe random and aperiodic metasurfaces [51] beyond the quasi-periodicity and phase discretization. Furthermore, ME is efficient both in memory and calculation time, allowing calculation of a large number of particles [52]. Therefore, ME can be used for optimization [53], inverse-design [54], and Monte-Carlo simulations. Metasurfaces typically consider planewave incidences, but interaction between metasurfaces and structured optical fields may result very different physics, such as recently investigated OAM-dependent metasurfaces [31].