Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group
Journal Home About Issues in Progress Current Issue All Issues Feature Issues

Efficient method for modeling large-scale arrays of optical nanoresonators based on the coupling theory of quasinormal mode

Qiyong Tao, Yuhang Su, Can Tao, Ying Zhong, and Haitao Liu

Open Access Open Access

Abstract

We propose an efficient method for calculating the electromagnetic field of a large-scale array of optical nanoresonators based on the coupling theory of quasinormal mode (QNM). In this method, two approaches of the scattered-field reconstruction and stationary-phase-principle calculated plane-wave expansion are developed to obtain the regularized QNM (RQNM) in different regions. This accurate and efficient calculation of RQNM resolves the far-field divergence issue of QNMs in the QNM-coupling theory, thus enabling a rapid computation of the electromagnetic field of a large-scale array of optical nanoresonators, which is a challenging task for full-wave numerical methods. Using this method, we consider the numerical example of the radiation problem of a single point source in a large-scale periodic array of optical nanoantennas. In comparison to full-wave numerical methods, this method significantly reduces the computation time by 1∼2 orders of magnitude while maintaining accuracy. The high computational efficiency and physical intuitiveness of the method enables to clarify the impact of array size (exceeding 50 × 50 wavelengths), period and field-coupling range (far beyond the tight-binding approximation) on the optical response. The proposed method and results can provide an efficient tool and guidance for the design of large-scale arrays of optical nanoresonators.

© 2024 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement

1. Introduction

Arrays of optical nanoresonators can simultaneously support localized resonance of individual resonators and mutual coupling between resonators, and thereby can achieve enhanced and flexibly controlled interactions between light and matter [1–5]. For instance, non-periodic arrays of optical nanoresonators can form metasurfaces [6–8] that enable flexible control of the amplitude, phase and polarization of the optical field, and have applications such as metalens [8] and holography [6]. Periodic arrays of optical nanoresonators support surface lattice resonance (SLR) [1,5,9] that enables large-area electromagnetic field enhancement [1,10,11], enhanced spontaneous emission [9] and far-field directional radiation [2,3], and have applications such as high-efficiency photovoltaic devices [4] and high-speed high-brightness light sources [1,12].

However, for solving the electromagnetic field of large-scale arrays of optical nanoresonators with an overall size much larger than the wavelength (reaching tens or even thousands of wavelengths), traditional full-wave numerical methods (such as finite-difference time-domain method [13] and finite element method [14]) suffer from enormous computational load which is even unachievable [15,16]. In recent years, quasinormal mode (QNM) expansion theory [17–20] in which the source-excited electromagnetic field is expressed as an expansion upon the basis of source-free eigensolutions of QNMs, shows high computational efficiency and physical intuitiveness, and has been widely applied in nanophotonics calculations [21], quantization of optical fields [22], design of optical nanoantennas [23,24], analysis of nonlinear effects in nanoresonators [25], and optical force calculations [26]. For coupled systems of optical nanoresonators, coupling theories of QNM also have significant developments in recent years [19,27–29]. However, difficulties arise when using the coupling theory of QNM to analyse large-scale arrays of optical nanoresonators due to the divergence of the far field of QNMs [19,27,29].

To address the aforementioned difficulties, this paper proposes an efficient method for calculating the electromagnetic field of large-scale arrays of optical nanoresonators based on the coupling theory of QNM. The QNM-coupling theory utilized in this paper [27] is rigorously derived from the first principle of Maxwell's equations, and therefore is generally applicable to coupled systems with energy loss and arbitrary forms of frequency dispersion of dielectric constant. In the near-field region, far-field region, and their inbetween transition region, we propose to use the QNM field, the plane-wave expansion calculated by stationary-phase principle, and the scattered-field reconstruction, respectively, to obtain the regularized QNM (RQNM). This accurate and efficient calculation of RQNM resolves the far-field divergence issue of QNMs in the QNM-coupling theory, thus enabling a rapid computation of electromagnetic field of large-scale arrays of optical nanoresonators—a task challenging for full-wave numerical methods. As numerical examples, we consider the radiation problem of a single point source in a large-scale periodic array of optical nanoantennas. In comparison to full-wave numerical methods, our method significantly reduces the computation time by 1∼2 orders of magnitude while maintaining accuracy. This high computational efficiency allows an efficient analysis of the impact of the array size and period on the optical response. Benefiting from the intuitive physical insight into the interaction between scatterers provided by the QNM-coupling theory, we quantitatively clarify the field-coupling range required to form the SLR in periodic array, which is found up to 31 times of the wavelength and far beyond the tight-binding approximation [30]. The proposed method and results can be used for an efficient design of large-scale arrays of optical nanoresonators [1–3,12].

2. Method

According to the coupling theory of QNM [27], for a coupled system consisting of P individual scatterers, the scattered field Ψs(r,ω) = [Es,Hs] can be expressed as a superposition of the scattered fields of each individual scatterer,

$${{\mathbf \Psi }^s}({\textbf r},\omega ) = \sum\limits_{p = 1}^P {{\mathbf \Psi }_p^s({\textbf r},\omega )} ,$$
where r = (x,y,z) represents spatial coordinates, and ω denotes the real-valued excitation angular frequency. In Eq. (1), the scattered field of the p-th scatterer is expanded upon the basis of its QNM field ${\tilde{\mathbf{\Psi}}_{p,m}}({\textbf r}) = [{\tilde{{\textbf E}}_{p,m}},{\tilde{{\textbf H}}_{p,m}}]$,
$${\mathbf \Psi }_p^s({\textbf r},\omega ) = \sum\limits_{m = 1}^{{M_p}} {{\alpha _{p,m}}(\omega ){{\tilde{\mathbf{\Psi}}}_{p,m}}({\textbf r})} .$$
In Eq. (2), the m-th QNM of the p-th scatterer is the eigensolution of the source-free Maxwell’s equtions at a complex eigenfrequency ${\tilde{\omega }_{p,m}}$,
$$ \nabla \times \tilde{\mathbf{E}}_{p, m}=i \tilde{\omega}_{p, m} \mu_0 \tilde{\mathbf{H}}_{p, m}, \nabla \times \tilde{\mathbf{H}}_{p, m}=-i \tilde{\omega}_{p, m} \boldsymbol{{\mathrm{\varepsilon}}}_p\left(\mathbf{r}, \tilde{\omega}_{p, m}\right) \cdot \tilde{\mathbf{E}}_{p, m}, $$
and satisfies the outgoing-wave condition at infinity. In Eq. (3), μ0 is the permeability in vacuum for non-magnetic medium, and ɛp is the permittivity tensor with the presence of only the p-th scatterer. After obtaining the QNMs by solving Eq. (3), the QNM expansion coefficients αp,m(ω) in Eq. (2) can be determined by solving the following set of linear equations:
$$[{\omega {\textbf I} - {\textbf K}(\omega )} ]{\textbf a}(\omega ) = \textbf{b}(\omega ).$$
In Eq. (4), the (p,m)-th row of the column vector a is αp,m, and I is the identity matrix. The element κ(p,m),(q,n) of the (p,m)-th row and (q,n)-th column of the matrix K is the QNM-coupling coefficient defined as (for p≠q):
$$ \kappa_{(p, m),(q, n)}(\omega)=\frac{-\omega}{F_{p, m}} \iiint_{V_p} \tilde{\mathbf{E}}_{p, m}(\mathbf{r}) \cdot \Delta \boldsymbol{{\mathrm{\varepsilon}}}_p(\mathbf{r}, \omega) \cdot \tilde{\mathbf{E}}_{q, n}(\mathbf{r}) d^3 \mathbf{r} $$
where Δɛp=ɛp−ɛb is the difference between the permittivity ɛp of the p-th scatterer and the permittivity ɛb of the background medium in absence of any scatterer. For p = q, κ(p,m),(q,n) is defined as ${\delta _{m,n}}{\tilde{\omega }_{p,m}}$. The integral region Vp is the region of the p-th scatterer, and Fp,m is the pseudoenergy of the QNM expressed as [18,31]:
$$ F_{p, m}=\iiint_{R^3}\left\{\left.\tilde{\mathbf{E}}_{p, m}(\mathbf{r}) \cdot \frac{\partial\left[\omega \boldsymbol{{\mathrm{\varepsilon}}}_p(\mathbf{r}, \omega)\right]}{\partial \omega}\right|_{\omega=\tilde{\omega}_{p, m}} \cdot \tilde{\mathbf{E}}_{p, m}(\mathbf{r})-\mu_0 \tilde{\mathbf{H}}_{p, m}(\mathbf{r}) \cdot \tilde{\mathbf{H}}_{p, m}(\mathbf{r})\right\} d^3 \mathbf{r}, $$
where R3 represents the entire space. In Eq. (4), the (p,m)-th row of the column vector b is the QNM excitation coefficient βp,m defined as:
$$ \beta_{p, m}(\omega)=\frac{-\omega}{F_{p, m}} \iiint_{V_p} \tilde{\mathbf{E}}_{p, m}(\mathbf{r}) \cdot \Delta \boldsymbol{{\mathrm{\varepsilon}}}_p(\mathbf{r}, \omega) \cdot \mathbf{E}^b(\mathbf{r}, \omega) d^3 \mathbf{r}, $$
where Eb is the background field as the excitation source, i.e., the incident field in the background medium in absence of any scatterer. Substituting the a solved from Eq. (4) into Eq. (2) and then into Eq. (1), one can obtain the scattered field of the coupled system. By superimposing the scattered field Ψs and the background field Ψb, the total field Ψ of the system is obtained, i.e., Ψ=Ψs+Ψb.

For large-scale arrays of scatterers, when calculating the QNM-coupling coefficients κ(p,m),(q,n) in Eq. (5), it is required to precalculate the QNM fields over a spatial range that extends to the whole size of the array. Therefore, when the array size is much larger than the wavelength λ, full-wave numerical calculations of QNM fields by solving Eq. (3) become computationally expensive and are even unachievable. A more serious issue is that, as |r|→∞, the QNM field $|{{{\tilde{\mathbf{\Psi}}}_{p,m}}({\textbf r})} |\to \infty$ [19,21,27,29], which leads to an exponential divergence of κ(p,m),(q,n) with the increase of the distance between scatterers.

To address the above issues in the QNM-coupling theory, we replace the QNM fields ${\tilde{\mathbf{\Psi}}_{p,m}}({\textbf r})$ in Eqs. (2) and (5) with the RQNM fields $\tilde{\mathbf{\Psi}}_{p,m}^R({\textbf r},\omega )$ [19,22] at a given real-valued excitation frequency ω. Within the near-field region (containing the scatterer and with a boundary S) where the QNM fields do not diverge, the RQNM fields are chosen to be the QNM fields. Outside S, the RQNM fields should satisfy the Maxwell's equations at the real frequency ω (hence depending on ω) and therefore have no divergence. On both sides of S, the RQNM fields should satisfy the continuity of tangential components of the electric and magnetic fields.

Firstly, outside the near-field region (i.e., outside the aforementioned S), we develop an approach to reconstruct the RQNM using the scattered fields from different excitation sources. According to Eq. (2), for the p-th scatterer we can set Mp excitation sources ${\textbf J}_p^{(\mu )}({\textbf r},\omega )$ (μ=1, 2, …, Mp) with the same real frequency ω but different distributions. Then the corresponding scattered fields ${\mathbf \Psi }_{p,\mu }^s({\textbf r},\omega )$ can be expanded with RQNMs:

$${\mathbf \Psi }_{p,\mu }^s({\textbf r},\omega ) = \sum\limits_{m = 1}^{{M_p}} {\alpha _{p,m}^{(\mu )}(\omega )\tilde{\mathbf{\Psi}}_{p,m}^R({\textbf r},\omega )} ,$$
where $\alpha _{p,m}^{(\mu )}(\omega )$ is the corresponding QNM expansion coefficient [18,31]:
$$ \alpha_{p, m}^{(\mu)}(\omega)=\left[-\omega \iiint_{V_p} \tilde{\mathbf{E}}_{p, m}(\mathbf{r}) \cdot \Delta\boldsymbol{{\mathrm{\varepsilon}}}_p(\mathbf{r}, \omega) \cdot \mathbf{E}_p^{(\mu)}(\mathbf{r}, \omega) d^3 \mathbf{r}\right] /\left[\left(\omega-\tilde{\omega}_{p, m}\right) F_{p, m}\right], $$
with ${\textbf E}_p^{(\mu )}({\textbf r},\omega )$ representing the background field excited by the ${\textbf J}_p^{(\mu )}({\textbf r},\omega )$ in absence of the scatterer. The scattered fields ${\mathbf \Psi }_{p,\mu }^s({\textbf r},\omega )$ can be obtained through full-wave numerical calculations. Thus, Eq. (7) forms a set of linear equations satisfied by the RQNM fields $\tilde{\mathbf{\Psi}}_{p,m}^R({\textbf r},\omega )$, whose solution then gives $\tilde{\mathbf{\Psi}}_{p,m}^R({\textbf r},\omega )$. The non-divergence of ${\mathbf \Psi }_{p,\mu }^s({\textbf r},\omega )$ at real ω ensures the non-divergence of the solved $\tilde{\mathbf{\Psi}}_{p,m}^R({\textbf r},\omega )$.

However, for the above approach, as the spatial range of the scattered fields obtained by full-wave computation increases, the computational amount increases rapidly and becomes even unachievable (for instance, when the spatial range reaches tens or even hundreds of wavelengths). In response to this issue, we develop a second approach for a rapid calculation of the RQNM in the far-field region. In this approach, utilizing the non-divergent QNM field on the boundary S of near-field region mentioned earlier, its plane-wave expansion integral in the far-field region at the given real frequency ω then gives the RQNM far field without divergence and can be calculated analytically using the stationary-phase principle [32,33]. The obtained expression for the m-th RQNM field of the p-th scatterer is (the derivation process can be found in Supplement 1 Sec. S1):

$$\tilde{\mathbf{\Psi}}_{p,m}^R({\textbf r},\omega ) = \sum\limits_{\sigma = 1,2} {\tilde{C}_{p,m}^ + ({{\tilde{k}}_x},{{\tilde{k}}_y},\sigma ){\mathbf \Psi }_{\textrm{PW}}^ + ({{\tilde{k}}_x},{{\tilde{k}}_y},\sigma )\textrm{exp}(i{k_0}{n_b}r)} \frac{\lambda }{{i{n_b}r}},$$
where the plane wave expansion coefficient is given by,
$$\tilde{C}_{p,m}^ + ({\tilde{k}_x},{\tilde{k}_y},\sigma ) = \frac{{k_0^2{n_b}{\eta _{\textrm{vac}}}\mathop{{\int\!\!\!\int}\mkern-11.4mu\circ}\nolimits_{S} {({\textbf E}_{\textrm{PW}}^ -{\times} {{\tilde{{\textbf H}}}_{p,m}} - {{\tilde{{\textbf E}}}_{p,m}} \times {\textbf H}_{\textrm{PW}}^ - )} \cdot \hat{{\textbf n}}ds }}{{8{\pi ^2}[{{\textbf E}_{\textrm{PW}}^ - ( - {{\tilde{k}}_x}, - {{\tilde{k}}_y},\sigma ) \cdot {\textbf E}_{\textrm{PW}}^ + ({{\tilde{k}}_x},{{\tilde{k}}_y},\sigma )} ]}}.$$
In Eqs. (8) and (9), ηvac is the wave impedance in vacuum, k0 = 2π/λ is the wave number (λ being the wavelength), nb is the refractive index of the uniform background medium, and $({\tilde{k}_x},{\tilde{k}_y},{\tilde{k}_z}) = {\textbf k}/({k_0}{n_b}) = (\sin \theta \cos \phi ,\sin \theta \sin \phi ,\cos \theta )$ are the direction cosines of the plane wave wavevector k, with $r \in [0, + \infty )$, $\theta \in [0,\pi ]$, $\phi \in [0,2\pi )$ representing the radial coordinate, polar angle, and azimuthal angle of r in spherical coordinates (coordinate origin at the center of the p-th scatterer), and (x,y,z) = (rsinθcosϕ,rsinθsinϕ,rcosθ) representing the Cartesian coordinates of r. ${\mathbf \Psi }_{\textrm{PW}}^ + ({\tilde{k}_x},{\tilde{k}_y},\sigma )\textrm{exp}[i{k_0}{n_b}({\tilde{k}_x}x + {\tilde{k}_y}y + {\tilde{k}_z}z)]$ denotes the up-going plane wave component, where σ=1, 2 corresponds to TE and TM polarizations, respectively, with the electric-field vector ${\textbf E}_{\textrm{PW}}^ +$ and magnetic-field vector ${\textbf H}_{\textrm{PW}}^ +$ perpendicular to the plane determined by k and z-axis, respectively. $[{\tilde{{\textbf E}}_{p,m}},{\tilde{{\textbf H}}_{p,m}}]$ represents the non-divergent QNM field on S. $[{\textbf E}_{\textrm{PW}}^ - ,{\textbf H}_{\textrm{PW}}^ - ]$ denotes the electric and magnetic fields of the down-going plane wave ${\mathbf \Psi }_{\textrm{PW}}^ - ( - {\tilde{k}_x}, - {\tilde{k}_y},\sigma )\textrm{exp}[ - i{k_0}{n_b}({\tilde{k}_x}x + {\tilde{k}_y}y + {\tilde{k}_z}z)]$ on S. The application of Eqs. (8) and (9) requires two conditions: i) k0nbr≫1, i.e., the distance r from the field point to the scatterer is much greater than the wavelength λ; ii) the scatterer is situated in a background medium with a uniform refractive index (nb). Equations (8) and (9) give the far field of the RQNM outside S in the z > 0 half-space at a real frequency ω, which is no longer divergent. The far field of RQNM in the z < 0 half-space can be obtained in a similar way.

In comparison to existing methods in the literature, such as using the surface integral [19] or volume integral [22] of Green's function to calculate the RQNM outside the scatterer, the above approach provides a fully-analytical expression for the far field of RQNM with a lower computational cost. Additionally, in Ref. [19], RQNM is applied to the QNM-coupling theory without excitation source (for solving the QNMs of the coupled system), whereas in this paper, RQNM is applied to the QNM-coupling theory with excitation source [27] (for solving the source-excited electromagnetic field of the coupled system).

To confirm that the RQNM calculated by the above approaches can accurately reproduce the scattered field of a single resonator, we consider a nanoantenna of a single gold nanowire in air (nb = 1) with a square-shaped cross-section (side length D = 40 nm) and a length L = 247 nm, as shown in Fig. 1(a). The full-wave numerical calculations of electromagnetic fields in this paper are performed using the Fourier modal method (FMM) [34–36] unless otherwise specified. The QNM is solved using the pole-search approach [37,38]. The frequency-dependent permittivity of gold is described by the Drude model $\varepsilon (\omega ) = {\varepsilon _0}{\varepsilon _\infty } - {\varepsilon _0}\omega _0^2/({\omega ^2} + i\omega \gamma )$ [39], where ɛ∞=8.842, 2πc/ω0 = 0.164 μm, 2πc/γ=20.689 μm are obtained through fitting experimental data from visible to far-infrared wavelengths [40], c is the speed of light in vacuum, and ɛ0 is the permittivity in vacuum. The scattered field of this nanowire antenna can be accurately described by a single QNM with a high quality factor, whose field distribution is shown in Fig. 1(b). Its complex eigenfrequency is ${\tilde{\omega }_0}/(2\pi c) = (\textrm{0}\textrm{.8842}-\mathrm{0.0477i)\ \mathrm{\mu} }{\textrm{m}^{ - 1}}$, and its quality factor is ${Q_0} = \textrm{Re} ({\tilde{\omega }_0})/[ - 2{\mathop{\rm Im}\nolimits} ({\tilde{\omega }_0})] = 9.27$. Figures 1(c1)-(c5) plot the scattered field of the antenna in z = 0 plane at different azimuthal angles ϕ=90°, 60°, 30°, 15°, 0°. The excitation field is an x-polarized plane wave propagating along negative z-direction satisfying the normalization of Ex = 1 at the top surface of the antenna (z = D/2). The results of full-wave FMM calculation and QNM expansion [Eq. (2)] are represented by the red solid lines and blue dash-dot lines, respectively. The black dashed lines are obtained by replacing the QNM field in the QNM expansion with the RQNM field obtained by the plane-wave expansion calculated by the stationary-phase principle (PESP). For the PESP, the S in Eq. (9) is selected to be a cuboid surface surrounding the antenna with a distance w = 50 nm to the antenna surface. As shown in Supplement 1 Fig. S1, the w that ensures the accuracy of PESP can be selected within a rather wide range. It is worth noting that, since the scattered field is only expanded with one QNM, this scattered field is simply the RQNM reconstructed by the scattered field. The results indicate that, in the near-field region (r=|r|<0.5λ), the QNM expansion is non-divergent with a high accuracy compared to the full-wave FMM calculation, while the PESP shows a large error. In the far-field region (r > 5λ), the PESP agrees well with the full-wave FMM calculation, while the QNM expansion exhibits significant divergence with a large error. It is noteworthy that for ϕ=0° [Fig. 1(c5)], the scattered far field obtained from the full-wave FMM calculation drops below 10−4 of the near field, and the field calculated with the PESP is theoretically 0 (numerically, |Ex| down to 1.9 × 10−19 at r = 0.01λ) due to the mirror symmetry of QNM electric field about y = 0 and z = 0, so that the curve of PESP becomes invisible. This implies that for ϕ=0°, although the PESP cannot precisely reproduce the very-weak scattered far field, the resultant error in the QNM-coupling theory can be neglected. In the transition region between the near-field and far-field regions (0.5λ<r < 5λ), the QNM expansion still exhibits divergence with a distinct error compared to the full-wave FMM calculation, and the PESP also shows some error for small values of ϕ (see ϕ=30°, 15°).

 figure: Fig. 1.

Fig. 1. (a) Schematic diagram of the antenna of a single gold nanowire with the origin O located at the center of the antenna. (b) Distribution of Re(Ex) (top) and Re(Ey) (bottom) in z=0 plane for the electric-field x and y components of the QNM which is normalized to have Ex(L/2+10 nm,0,0)=1. (c1)-(c5) Distribution of |Ex| in z=0 plane at azimuthal angles ϕ=90°, 60°, 30°, 15°, and 0° for the scattered field of the antenna under the excitation of an x-polarized plane wave propagating along the negative z-direction. In the horizontal axis, r=|r| and the wavelength is λ=1127 nm. The red solid line, blue dash-dot line, and black dashed line represent the results obtained by full-wave FMM calculation, QNM expansion, and PESP, respectively. The insets show local enlargements within the range of r/λ<0.3.

Download Full Size | PDF

Based on the above results, in the subsequent numerical implementation of the QNM-coupling theory, the RQNM field in the near-field region (r < 0.5λ) is selected as the non-divergent QNM field. In the transition region (0.5λ<r < 5λ), the RQNM field is obtained from the scattered-field reconstruction. In the far-field region (r > 5λ), the RQNM field is obtained with the PESP approach.

For the calculation of the electromagnetic field of large-scale arrays of optical nanoresonators, compared to computationally demanding full-wave numerical methods, the proposed QNM-coupling method has the following advantages in reducing computational load:

  • i) For the calculation of RQNM of individual resonators, the sizes of near-field and transition regions are on the order of the wavelength and much smaller than the size of the entire array. So the computational load of the full-wave calculation of RQNM in the near-field and transition regions is feasible and much lower than that for the entire array. In the far-field region whose size is much larger than the wavelength, the RQNM is expressed analytically by the stationary-phase-principle calculated plane-wave expansion which requires a small computational effort.
  • ii) The calculation of the coefficient matrix K of the linear equations (4) for solving the QNM expansion coefficients and the solution of the linear equations require significantly less computation compared to full-wave electromagnetic field calculations of the array.
  • iii) When changing the number or relative positions of resonators in the array (for example, altering the period of the periodic array or changing the arrangement of resonators in a non-periodic array), as long as the individual resonators remain unchanged, there is no need to recalculate the RQNM for calculating the electromagnetic field of the array.
  • iv) When changing the excitation source of background field Eb (for example, altering its polarization state or spatial distribution, but keeping its real-valued excitation frequency fixed), there is no need to recalculate the RQNM and its coupling coefficients κ(p,m),(q,n) for calculating the electromagnetic field of the array.

3. Result

To demonstrate the advantages of our method in reducing computational amount and providing physical insight, we consider large-scale periodic arrays of nanoresonators [2,3,9] as numerical examples. It is important to note that our method is applicable as well to non-periodic arrays of resonators, such as metalens [8] and metasurface holographic devices [6]. Leveraging the symmetry of the periodic array of resonators, the computational cost of the RQNM-coupling coefficients κ(p,m),(q,n) in Eq. (5) can be further reduced. Assuming that the scattered field of each scatterer of resonator is expanded by M RQNMs, for a periodic array composed of Px rows and Py columns of identical scatterers (total number of scatterers P = Px × Py), the total number of κ(p,m),(q,n) is (PxPyM)2. For two scatterers with unchanged relative positions (but with possibly changed absolute positions in the array), their κ(p,m),(q,n) should be unchanged, which reduces the number of different-valued κ(p,m),(q,n) to (2Px−1)(2Py−1)M2. If the RQNM fields exhibit symmetry in both x and y directions (such as the single-nanowire antenna and dipole nanoantenna in the numerical examples), the number of different-valued κ(p,m),(q,n) is further reduced to PxPyM2.

To numerically test the accuracy of the proposed QNM-coupling method, we consider a periodic array of nanoantennas under excitation by a single point source [9,12,41], and compare the results obtained from full-wave FMM numerical calculations and those from the QNM-coupling method. The parameters of an individual antenna are the same as those in Fig. 1, and the array periods in x and y directions are all 1125 nm, as shown in Fig. 2(a). The x-polarized electric-current point source can be expressed as the electric-current density ${\textbf J}({\textbf r}) = \hat{{\textbf x}}\delta ({\textbf r} - {{\textbf r}_0})$, with a radiation wavelength of λ=1127 nm. The source is located to the right of the central antenna and 20 nm away from its end face, i.e., at r0 = (x0,y0,z0) = (L/2 + 20 nm,0,0), where L = 247 nm is the antenna length, and the coordinate origin O is set at the center of the central antenna. As illustrated in Fig. 1, the optical response of each antenna in the array is dominated by only one QNM in the QNM expansion, and therefore only this dominant QNM is considered in the QNM-coupling method.

 figure: Fig. 2.

Fig. 2. (a) Schematic diagram of a periodic array of single-nanowire antennas (array period a=b=1125 nm), where the red dot represents the point source (emission wavelength λ=1127 nm). (b) For an array size of 5×5, the far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$ calculated by the full-wave FMM (b1) and the QNM-coupling method (b2). (c) Positions of Rayleigh anomalies in the plane of $({\tilde{k}_x},{\tilde{k}_y})$, where the values of (m,n) in Eq. (13) are also provided. (d1)-(d3) $P({\tilde{k}_x},{\tilde{k}_y})$ calculated by the QNM-coupling method for array sizes of 15×15, 31×31 and 51×51, respectively. (e) For an infinite periodic array, $P({\tilde{k}_x},{\tilde{k}_y})$ calculated by the reciprocity-theorem method. (f1) Far-field emission-rate enhancement factor γrad=Γrad/Γvac calculated by the QNM-coupling method for different numbers of array periods (blue solid line), and the γrad for an infinite periodic array calculated by the reciprocity-theorem method (red dashed line). (f2) Total emission-rate enhancement factor γtotal=Γtotal/Γvac calculated by the QNM-coupling method for different numbers of array periods (blue solid line), and the γtotal for an infinite periodic array calculated by the ASM (red dashed line). The numbers of array periods along the x and y directions are set to be the same. (g) Distribution of the electric-field x-component Re(Ex) (normalized by Γvac) around the central antenna in the z=0 plane. The results are obtained by the QNM-coupling method for an array size of 51×51 (g1) and by the ASM for an infinite periodic array (g2).

Download Full Size | PDF

For the full-wave FMM calculations, in the region of z > 0 (upper half-space), the far-field time-averaged Poynting vector ${\textbf S}({\tilde{k}_x},{\tilde{k}_y},r)$ is given by:

$${\textbf S}({\tilde{k}_x},{\tilde{k}_y},r) = {\left( {\frac{\lambda }{{{n_b}r}}} \right)^2}\sum\limits_{\sigma = 1,2} {|{{\tilde{C}}^ + }({{\tilde{k}}_x},{{\tilde{k}}_y},\sigma ){|^2}{\textbf S}_{\textrm{PW}}^ + ({{\tilde{k}}_x},{{\tilde{k}}_y},\sigma )} ,$$
$${\textbf S}_{\textrm{PW}}^ + ({\tilde{k}_x},{\tilde{k}_y},\sigma ) = \frac{1}{2}\textrm{Re} [{\textbf E}_{\textrm{PW}}^ + ({\tilde{k}_x},{\tilde{k}_y},\sigma ) \times {\textbf H}_{\textrm{PW}}^ + {({\tilde{k}_x},{\tilde{k}_y},\sigma )^\ast }],$$
where the up-going plane-wave expansion coefficient ${\tilde{C}^ + }({\tilde{k}_x},{\tilde{k}_y},\sigma )$ can be obtained by replacing the ${\tilde{\mathbf{\Psi}}_{p,m}} = [{\tilde{{\textbf E}}_{p,m}},{\tilde{{\textbf H}}_{p,m}}]$ in Eq. (9) with the total field Ψ=[E,H] of the periodic array excited by a single point source. Here Ψ is obtained from FMM calculations, and the closed surface S is chosen to enclose the entire antenna array (with a distance of 300 nm to the boundary of the array). ${\textbf S}_{\textrm{PW}}^ + ({\tilde{k}_x},{\tilde{k}_y},\sigma )$ is the Poynting vector of the plane-wave mode with unitary coefficient. An enhancement factor of the far-field radiation intensity is defined as $P({\tilde{k}_x},{\tilde{k}_y}) = {{|{\textbf S}({{\tilde{k}}_x},{{\tilde{k}}_y},r)|} / {{S_{\textrm{vac}}}}}$, where Svac=Γvac/(4πr2) is an average radiation intensity of a point source in vacuum on a sphere of radius r, and $\Gamma_{\mathrm{vac}}=\eta_{\mathrm{vac}} k_0^2 n_{\mathrm{vac}} /(12 \pi)$ is the radiation power of the point source in vacuum, with ηvac and nvac = 1 being the wave impedance and refractive index in vacuum, respectively. It is notable that $P({\tilde{k}_x},{\tilde{k}_y})$ is independent of r and is therefore also called the far-field radiation-intensity angular distribution.

Using the QNM-coupling method, one can obtain for Eq. (10) (proof in Supplement 1 Sec. S2):

$${\tilde{C}^ + }({\tilde{k}_x},{\tilde{k}_y},\sigma ) = {\tilde{C}^b}({\tilde{k}_x},{\tilde{k}_y},\sigma ) + \sum\limits_{p = 1}^P {\sum\limits_{m = 1}^{{M_p}} {{\alpha _{p,m}}\tilde{C}_{p,m}^ + ({{\tilde{k}}_x},{{\tilde{k}}_y},\sigma )g({{\tilde{k}}_x},{{\tilde{k}}_y},{{\textbf r}_p})} } ,$$
where ${\tilde{C}^b}({\tilde{k}_x},{\tilde{k}_y},\sigma )$ is the up-going plane-wave expansion coefficient for the background field Ψb excited by a point source in the homogeneous free space [expression given in Supplement 1 Eq. (S2.13)], $\tilde{C}_{p,m}^ + ({\tilde{k}_x},{\tilde{k}_y},\sigma )$ is given by Eq. (9), $g({\tilde{k}_x},{\tilde{k}_y},{{\textbf r}_p}) = \textrm{exp} ( - i{k_0}{n_b}{r_p}\cos {\alpha _p})$, rp is the position vector of the center ${O^{\prime}_p}$ of the p-th scatterer, αp is the angle between the field position vector r and rp [as shown in Fig. 2(a)], rp=|rp|, and the coordinate origin O is set at the center of the central scatterer. In spherical coordinates, r = (rsinθcosϕ,rsinθsinϕ,rcosθ), rp = (rpcosϕp,rpsinϕp,0), and therefore cosαp = sinθcos(ϕ−ϕp). Substituting Eq. (12) into Eq. (10) then allows the calculation of the far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$.

First, we consider a 5 × 5 array of the single-nanowire antennas. Figures 2(b1)-(b2) show the $P({\tilde{k}_x},{\tilde{k}_y})$ in the z > 0 region obtained by both full-wave FMM numerical calculations and the QNM-coupling method. Their good agreement verifies the validity of the QNM-coupling method.

Next, to further confirm the validity of the QNM-coupling method, we continue to increase the array size. It is expected that when the array size is large enough, the field distribution of the source radiation will approach the distribution in an infinite periodic array. For the calculation of the radiation field of a single point source in an infinite periodic array, the periodicity of the electromagnetic field is broken by the single point source, making it unfeasible to calculate the field directly. In this case, the reciprocity-theorem method [41–44] can be employed to calculate the far-field distribution (details of this method can be found in Supplement 1 Sec. S3), and the array scanning method (ASM) [41,45] can be used to calculate the near-field distribution (details of this method can be found in Supplement 1 Sec. S4).

Figure 2(d) presents the far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$ calculated with the QNM-coupling method for gradually increased array sizes. The sampling intervals for polar angle θ and azimuthal angle ϕ are both 1°. Figure 2(e) shows the $P({\tilde{k}_x},{\tilde{k}_y})$ obtained by the reciprocity-theorem method for an infinite periodic array of antennas excited by a point source. The calculation of Fig. 2(e) is performed with the full-wave FMM, and due to computational-amount constraints, the sampling intervals for θ and ϕ are both 2°. Comparison between Figs. 2(d) and 2(e) shows that when the array size reaches 15 × 15, the $P({\tilde{k}_x},{\tilde{k}_y})$ is already similar to the result of the infinite periodic array. When the array size reaches 51 × 51, the $P({\tilde{k}_x},{\tilde{k}_y})$ is nearly identical to the result of the infinite periodic array, thereby confirming the validity of the QNM-coupling method.

The positions of Rayleigh anomaly (RA) in the plane of $({\tilde{k}_x},{\tilde{k}_y})$ are shown in Fig. 2(c) and are given by:

$${({k_x} - m2\pi /a)^2} + {({k_y} - n2\pi /b)^2} = {({k_0}{n_b})^2},$$
where $({k_x},{k_y}) = {k_0}{n_b}({\tilde{k}_x},{\tilde{k}_y})$ are the x and y components of the wave vector of the far-field radiation plane wave, nb = 1 is the refractive index of the air free space above the antenna, k0 is the wave number in vacuum, a and b are the array periods in x and y directions, respectively (a = b = 1125 nm in Fig. 2). The results in Fig. 2(d)-(e) indicate that there is a sudden change in $P({\tilde{k}_x},{\tilde{k}_y})$ near the RA [9,46], and this sudden change becomes more pronounced for larger array size.

Figure 2(f1)-(f2) present the curves of the far-field emission-rate enhancement factor γrad=Γrad/Γvac and the total emission-rate enhancement factor (i.e., Purcell factor) γtotal=Γtotal/Γvac of the point source as a function of the number of periods in x and y directions. Here Γvac is the radiation power of the point source in vacuum [expression given after Eq. (11)]. The far-field emission rate is expressed as ${\Gamma _{\textrm{rad}}} = \mathop{{\int\!\!\!\int}\mkern-11.4mu\circ}\nolimits_A {\hat{{\textbf n}} \cdot {\textbf S}da}$, where A is a closed surface encompassing the entire array, $\hat{{\textbf n}}$ is the outward normal unit vector on A, and S is the time-averaged Poynting vector. The x-polarized electric-current point source can be expressed as the electric-current density ${\textbf J}({\textbf r}) = \hat{{\textbf x}}\delta ({\textbf r} - {{\textbf r}_0})$, whose total spontaneous emission rate can be expressed as Γtotal = −Re[Ex(r0)]/2, with Re[Ex(r0)] representing the real part of the x-component of the electric field at the position of the point source. Increasing γrad can enhance fluorescence intensity, which is useful for single-molecule fluorescence detection [47] and high-brightness light sources [3]. Elevating γtotal can improve the modulation rate of light source, which is beneficial for high-speed communication light sources such as light-emitting diodes [48] and single-photon sources [49]. The results in Fig. 2(f) indicate that when the number of array periods exceeds 10, the values of γrad and γtotal computed by the QNM-coupling method stabilize and agree with the results of an infinite periodic array obtained using the reciprocity-theorem method and ASM (both employing the full-wave FMM). As shown in Fig. 2(g), the QNM-coupling method agrees well with the ASM in calculating the near field around the central antenna.

In the above calculations, the excitation frequency ω is close to the resonant frequency [i.e., the eigenfrequency’s real part $\textrm{Re} ({\tilde{\omega }_0})$ of the dominant QNM] of each antenna in the array. As shown in Supplement 1 Fig. S2, when ω is moderately away from $\textrm{Re} ({\tilde{\omega }_0})$ (at ω1 or ω2 in Fig. S2), the QNM-coupling method considering the dominant QNM can still accurately calculate the $P({\tilde{k}_x},{\tilde{k}_y})$. When ω is significantly away from $\textrm{Re} ({\tilde{\omega }_0})$, the QNM-coupling method is expected to provide accurate predictions but at the expense of considering more QNMs to describe the optical response of each resonator in the array [18]. These results demonstrate that for large-scale periodic arrays of nanoresonators, the QNM-coupling method can provide highly accurate results for both the near and the far fields.

Compared with full-wave numerical methods, the QNM-coupling method exhibits a significant advantage in reducing computational amount [27,50]. The computational times for the main results in this study are presented in Table 1. The computational time of the QNM-coupling method is mainly devoted to solving the QNM and RQNM fields. Once these mode fields are obtained, the electromagnetic field for a large-scale array, such as the 51 × 51 array calculated in Fig. 2(d3) with a total size exceeding 50 × 50 wavelengths, can be obtained in just 73 seconds. Such a computation would be challenging if using full-wave numerical methods. In comparison with full-wave numerical methods, the QNM-coupling method reduces the computational time by 1 to 2 orders of magnitude.

Tables Icon

Table 1. Computation times of different methods. The computation is performed on a desktop PC with 6-core CPU of 4 GHz and a 64 GB RAM. The total number of sampling points of (θ,ϕ) for the reciprocity-theorem method in Fig. 2(e) and Fig. 4(d) is 1/4 of that for the QNM-coupling method. Therefore, the computation time given in the table for the reciprocity-theorem method is four times of the actual computation time.

View Table

Benefiting from the high computational efficiency, the QNM-coupling method enables an efficient design of large-scale arrays of nanoresonators. For example, as shown in Supplement 1 Fig. S3, by repeatedly using the already calculated QNM and RQNM fields, it takes only 9 min 44 s to obtain the far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$ for eight different array periods (array size 51 × 51). The results indicate that when choosing an array period of 1 μm, $P({\tilde{k}_x},{\tilde{k}_y})$ reaches a maximum value of approximately 1100 near θ=0. This result can be useful for enhancing the fluorescence collection efficiency of objectives, which is significant for applications such as high-brightness light sources [1,3] and highly sensitive fluorescence molecular sensing [51].

In addition, compared to the full-wave numerical method, the QNM-coupling method can provide an intuitive physical insight into the mutual coupling between resonators in the array. For example, Eq. (5) indicates that when the distance between two scatterers is sufficiently large, the QNM-coupling coefficient κ(p,m),(q,n) approaches 0 and thus can be neglected. For a quantitative analysis, in the following calculations, when the distance between the centers of any two scatterers is greater than RK, the corresponding QNM-coupling coefficient is artificially set to be κ(p,m),(q,n) = 0, where RK is called designated coupling distance. By changing the value of RK, a quantitative analysis of the influence of RK on the optical properties of the resonator array can be performed. This analysis is unachievable if using the full-wave numerical method. Figure 3 presents the far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$ for a 51 × 51 periodic array of single-nanowire antennas [with parameters as specified in Fig. 2(d3)] for different values of RK. The results show that when RK = 9λ, the positions where $P({\tilde{k}_x},{\tilde{k}_y})$ changes suddenly can roughly reflect the locations of Rayleigh anomaly [Eq. (13)]; when RK = 31λ, the result is almost identical to that of an infinite periodic array [Fig. 2(e)]. This indicates that to form SLR [5,9], the involved field-coupling distance between scatterers reaches 31 times of the wavelength. Therefore, the tight-binding approximation [30], which only considers the near-field coupling between adjacent scatterers, cannot reflect SLR. This conclusion holds as well when the array period is distinctly smaller or larger than the wavelength, as shown in Supplement 1 Fig. S4 (the period a = b = 1125 nm is close to the wavelength λ=1127 nm in Fig. 3).

 figure: Fig. 3.

Fig. 3. Far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$ for different values of the designated coupling distance RK. (a)-(h) RK=0 (i.e., setting QNM-coupling coefficients κ(p,m),(q,n)=0 for any p≠q), 5λ, 9λ, 15λ, 21λ, 31λ, 41λ, 51λ. The calculation is performed for a 51×51 periodic array of single-nanowire antennas, with parameters identical to those in Fig. 2(d3).

Download Full Size | PDF

In the calculations above, only one QNM with a dominant contribution was considered for each resonator in the array. For the more general case that multiple QNMs should be considered for each resonator in the array, the numerical test of the QNM-coupling method will be provided in the following. For this purpose, we consider the calculation of the radiation field of a single point source in a periodic array of dipole nanoantennas. The dipole nanoantenna is composed of two identical single-gold-nanowire antennas (parameters following those in Fig. 1) separated by an air gap of size d = 100 nm. As shown in Fig. 4(a), an x-polarized electric-current point source ${\textbf J}({\textbf r}) = \hat{{\textbf x}}\delta ({\textbf r} - {{\textbf r}_0})$ (red dot) is located to the left of the central antenna at a distance of 20 nm to the antenna end, i.e., at r0 = (−d/2 − L−20 nm,0,0), where L = 247 nm is the antenna arm length. The coordinate origin O is located at the center of the gap of the central antenna. The scattered field of this dipole nanoantenna should be described by two dominant QNMs with complex eigenfrequencies ${\tilde{\omega }_1}/(2\pi c) = (\textrm{0}\textrm{.9142} - \textrm{0}\mathrm{.0274i)\ \mathrm{\mu} }{\textrm{m}^{ - 1}}$ and ${\tilde{\omega }_2}/(2\pi c) = (\textrm{0}\textrm{.8486} - \textrm{0}\mathrm{.0668i)\ \mathrm{\mu} }{\textrm{m}^{ - 1}}$, corresponding to symmetric and antisymmetric modes, respectively. Their electric field distributions are shown in the upper and lower parts of Fig. 4(b). When solving RQNMs using Eq. (7), the two excitation sources are set as follows: 1) an x-polarized electric-current point source located at the center of the gap of the dipole antenna, i.e., r0 = (x0,y0,z0) = (0,0,0); 2) an x-polarized electric-current point source located to the left of the dipole antenna at a distance of 20 nm to the antenna end, i.e., r0 = (−d/2 − L−20 nm,0,0).

 figure: Fig. 4.

Fig. 4. (a) Schematic diagram of a periodic array of dipole nanoantennas, where the array period is a=b=1171 nm, and the point source (red dot) radiates at a wavelength λ=1171 nm. (b) Distributions in the z=0 plane for the electric field x-component Re(Ex) of symmetric QNM (upper image) and antisymmetric QNM (lower image) which are normalized to have Ex(d/2+L+10nm,0,0)=1. (c) Far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$ for a 51×51 periodic array of antennas calculated with the QNM-coupling method. (d) $P({\tilde{k}_x},{\tilde{k}_y})$ for an infinite periodic array of antennas calculated with the reciprocity-theorem method using the full-wave FEM.

Download Full Size | PDF

Figure 4(c) shows the far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$ calculated with the QNM-coupling method for a 51 × 51 periodic array of dipole antennas. Figure 4(d) provides the result of an infinite periodic array of dipole antennas calculated with the reciprocity-theorem method. Unlike the previous numerical results, significant differences are observed in the left and right directions, which reflects the contributions of multiple QNMs. The agreement between the QNM-coupling method and the reciprocity-theorem method confirms the validity of the QNM-coupling method for the general case of considering multiple QNMs for a single scatterer. The total computation time for the QNM-coupling method is 7.31 hours, which is much less than the 191 hours required for the reciprocity-theorem method (see Table 1). The reciprocity-theorem method is performed with the finite element method (FEM) which is executed with COMSOL Multiphysics software. The number of (θ,ϕ) sampling points for the reciprocity-theorem method is only 1/4 of that for the QNM-coupling method (the sampling intervals of θ and ϕ are all 2° for the former, and are all 1° for the latter).

4. Conclusion

We propose an efficient method for calculating the electromagnetic field of large-scale arrays of nanoresonators based on the coupling theory of QNM. In this method, we developed two approaches, namely, the scattered-field reconstruction and stationary-phase-principle calculated plane-wave expansion, to calculate the RQNM in the transition and far-field regions. This accurate and efficient calculation of RQNM resolves the far-field divergence issue of QNMs in the QNM-coupling theory, and thus enables an efficient modeling of large-scale arrays of optical nanoresonators—a task challenging for full-wave numerical methods. As numerical examples, we applied this method to calculate the radiation field of a single point source in a large-scale periodic array of nanantennas. Comparison with full-wave numerical methods indicates that our method can achieve a significant reduction in computation time by 1 to 2 orders of magnitude while maintaining accuracy. The high computational efficiency of the method enables to unveil the impact of the array size (exceeding 50 × 50 wavelengths) and period on the radiation properties. Utilizing the intuitive physical insights into the mutual coupling between nanoresonators provided by our method, we quantitatively determined the field-coupling range (RK ≥ 31λ) required for forming the SLR in a periodic array of nanoresonators. We emphasize that all these results are difficult to obtain if using full-wave numerical methods. Although the numerical examples in this paper consider scatterers in a homogeneous free space, the method remains applicable for scenarios with a substrate by further developing efficient approaches for calculating the far field of RQNM [32]. This method can serve as an efficient and intuitive tool for analysing nanophotonic devices composed of large-scale arrays of optical nanoresonators [1,6,8,9,12].

Funding

National Natural Science Foundation of China (62075104, 61775105).

Disclosures

The authors declare no conflicts of interest.

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

Supplemental document

See Supplement 1 for supporting content.

References

1. G. Lozano, S. R. Rodriguez, M. A. Verschuuren, et al., “Metallic nanostructures for efficient LED lighting,” Light: Sci. Appl. 5(6), e16080 (2016). [CrossRef]  

2. A. Vaskin, J. Bohn, K. E. Chong, et al., “Directional and spectral shaping of light emission with Mie-resonant silicon nanoantenna arrays,” ACS Photonics 5(4), 1359–1364 (2018). [CrossRef]  

3. S. Liu, A. Vaskin, S. Addamane, et al., “Light-emitting metasurfaces: simultaneous control of spontaneous emission and far-field radiation,” Nano Lett. 18(11), 6906–6914 (2018). [CrossRef]  

4. M. Van Lare, F. Lenzmann, M. A. Verschuuren, et al., “Mode coupling by plasmonic surface scatterers in thin-film silicon solar cells,” Appl. Phys. Lett. 101(22), 221110 (2012). [CrossRef]  

5. V. G. Kravets, A. V. Kabashin, W. L. Barnes, et al., “Plasmonic surface lattice resonances: a review of properties and applications,” Chem. Rev. 118(12), 5912–5951 (2018). [CrossRef]  

6. S. Zhang, H. Zhou, B. Liu, et al., “Recent advances and prospects of optical metasurfaces,” ACS Photonics 10(7), 2045–2063 (2023). [CrossRef]  

7. J. Hu, S. Bandyopadhyay, Y. Liu, et al., “A Review on metasurface: from principle to smart metadevices,” Front. Phys. 8, 586087 (2021). [CrossRef]  

8. S. So, J. Mun, J. Park, et al., “Revisiting the design strategies for metasurfaces: fundamental physics, optimization, and beyond,” Adv. Mater. 35(43), 2206399 (2023). [CrossRef]  

9. F. Laux, N. Bonod, and D. Gérard, “Single emitter fluorescence enhancement with surface lattice resonances,” J. Phys. Chem. C 121(24), 13280–13289 (2017). [CrossRef]  

10. G. Lozano, D. J. Louwers, S. R. Rodríguez, et al., “Plasmonics for solid-state lighting: enhanced excitation and directional emission of highly efficient light sources,” Light: Sci. Appl. 2(5), e66 (2013). [CrossRef]  

11. X. Zhang, H. Liu, and Y. Zhong, “Large-area electromagnetic enhancement by a resonant excitation of surface waves on a metallic surface with periodic subwavelength patterns,” Opt. Express 21(20), 24139 (2013). [CrossRef]  

12. A. Vaskin, R. Kolkowski, A. F. Koenderink, et al., “Light-emitting metasurfaces,” Nanophotonics 8(7), 1151–1198 (2019). [CrossRef]  

13. F. I. Baida and A. Belkhir, “Finite Difference Time Domain Method for Grating Structures,” in Gratings: Theory and Numeric Applications, E. Popov, ed., 2nd rev. ed. (Institut Fresnel, AMU2014), Chap. 9.

14. G. Demésy, F. Zolla, A. Nicolet, et al., eds., 2nd rev. ed. (Institut Fresnel, AMU2014), Chap. 5.

15. J. Skarda, R. Trivedi, L. Su, et al., “Low-overhead distribution strategy for simulation and optimization of large-area metasurfaces,” npj Comput. Mater. 8(1), 78 (2022). [CrossRef]  

16. T. W. Hughes, M. Minkov, V. Liu, et al., “A perspective on the pathway toward full wave simulation of large area metalenses,” Appl. Phys. Lett. 119(15), 150502 (2021). [CrossRef]  

17. E. S. C. Ching, P. T. Leung, A. Maassen Van Den Brink, et al., “Quasinormal-mode expansion for waves in open systems,” Rev. Mod. Phys. 70(4), 1545–1554 (1998). [CrossRef]  

18. P. Lalanne, W. Yan, K. Vynck, et al., “Light interaction with photonic and plasmonic resonances,” Laser Photonics Rev. 12(5), 1700113 (2018). [CrossRef]  

19. P. T. Kristensen, K. Herrmann, F. Intravaia, et al., “Modeling electromagnetic resonators using quasinormal modes,” Adv. Opt. Photonics 12(3), 612–708 (2020). [CrossRef]  

20. S. Both and T. Weiss, “Resonant states and their role in nanophotonics,” Semicond. Sci. Technol. 37(1), 013002 (2022). [CrossRef]  

21. C. Sauvan, T. Wu, R. Zarouf, et al., “Normalization, orthogonality, and completeness of quasinormal modes of open systems: the case of electromagnetism [Invited],” Opt. Express 30(5), 6846–6885 (2022). [CrossRef]  

22. S. Franke, S. Hughes, M. K. Dezfouli, et al., “Quantization of quasinormal modes for open cavities and plasmonic cavity quantum electrodynamics,” Phys. Rev. Lett. 122(21), 213901 (2019). [CrossRef]  

23. B. Gurlek, V. Sandoghdar, and D. Martín-Cano, “Manipulation of quenching in nanoantenna–emitter systems enabled by external detuned cavities: a path to enhance strong-coupling,” ACS Photonics 5(2), 456–461 (2018). [CrossRef]  

24. T. Wu, W. Yan, and P. Lalanne, “Bright plasmons with cubic nanometer mode volumes through mode hybridization,” ACS Photonics 8(1), 307–314 (2021). [CrossRef]  

25. C. Gigli, T. Wu, G. Marino, et al., “Quasinormal-mode non-Hermitian modeling and design in nonlinear nano-optics,” ACS Photonics 7(5), 1197–1205 (2020). [CrossRef]  

26. Z. Qi, C. Tao, S. Rong, et al., “Efficient method for the calculation of the optical force of a single nanoparticle based on the quasinormal mode expansion,” Opt. Lett. 46(11), 2658 (2021). [CrossRef]  

27. C. Tao, J. Zhu, Y. Zhong, et al., “Coupling theory of quasinormal modes for lossy and dispersive plasmonic nanoresonators,” Phys. Rev. B 102(4), 045430 (2020). [CrossRef]  

28. B. Vial and Y. Hao, “A coupling model for quasi-normal modes of photonic resonators,” J. Opt. 18(11), 115004 (2016). [CrossRef]  

29. K. G. Cognée, Hybridization of Open Photonic Resonators, Ph.D. dissertation (Université de Bordeaux; Universiteit van Amsterdam, 2020).

30. B. Xi, H. Xu, S. Xiao, et al., “Theory of coupling in dispersive photonic systems,” Phys. Rev. B 83(16), 165115 (2011). [CrossRef]  

31. C. Sauvan, J. P. Hugonin, I. S. Maksymov, et al., “Theory of the spontaneous optical emission of nanosize photonic and plasmon resonators,” Phys. Rev. Lett. 110(23), 237401 (2013). [CrossRef]  

32. J. Yang, J.-P. Hugonin, and P. Lalanne, “Near-to-far field transformations for radiative and guided waves,” ACS Photonics 3(3), 395–402 (2016). [CrossRef]  

33. N. Wang, Y. Zhong, and H. Liu, “Spontaneous emission enhancement by rotationally-symmetric optical nanoantennas: impact of radially and axially propagating surface plasmon polaritons,” Opt. Express 30(8), 12797–12822 (2022). [CrossRef]  

34. L. Li, “New formulation of the Fourier modal method for crossed surface-relief gratings,” J. Opt. Soc. Am. A 14(10), 2758–2767 (1997). [CrossRef]  

35. J. P. Hugonin and P. Lalanne, “Perfectly matched layers as nonlinear coordinate transforms: a generalized formalization,” J. Opt. Soc. Am. A 22(9), 1844–1849 (2005). [CrossRef]  

36. H. Liu, DIF CODE for Modeling Light Diffraction in Nanostructures (Nankai University, Tianjin, 2010).

37. Y. Li, H. Liu, H. Jia, et al., “Fully vectorial modeling of cylindrical microresonators with aperiodic Fourier modal method,” J. Opt. Soc. Am. A 31(11), 2459–2466 (2014). [CrossRef]  

38. P. Lalanne, W. Yan, A. Gras, et al., “Quasinormal mode solvers for resonators with dispersive materials,” J. Opt. Soc. Am. A 36(4), 686–704 (2019). [CrossRef]  

39. A. Raman and S. Fan, “Photonic band structure of dispersive metamaterials formulated as a Hermitian eigenvalue problem,” Phys. Rev. Lett. 104(8), 087401 (2010). [CrossRef]  

40. E. D. Palik, Handbook of Optical Constants of Solids (Academic press, New York, 1998).

41. Y. H. Su, L. Zhang, C. Tao, et al., “Spontaneous emission enhancement and directional emission by an optical nanonatenna array on a metallic mirror,” Acta Phys. Sin. 72(7), 078101 (2023). [CrossRef]  

42. O. T. A. Janssen, A. J. H. Wachters, and H. P. Urbach, “Efficient optimization method for the light extraction from periodically modulated LEDs using reciprocity,” Opt. Express 18(24), 24522–24535 (2010). [CrossRef]  

43. M. Kahl and E. Voges, “Analysis of plasmon resonance and surface-enhanced Raman scattering on periodic silver structures,” Phys. Rev. B 61(20), 14078–14088 (2000). [CrossRef]  

44. S. Zhang, E. R. Martins, A. G. Diyaf, et al., “Calculation of the emission power distribution of microstructured OLEDs using the reciprocity theorem,” Synth. Met. 205, 127–133 (2015). [CrossRef]  

45. F. Capolino, D. R. Jackson, D. R. Wilton, et al., “Comparison of methods for calculating the field excited by a dipole near a 2-D periodic material,” IEEE Trans. Antennas Propag. 55(6), 1644–1655 (2007). [CrossRef]  

46. G. Vecchi, V. Giannini, and J. G. Rivas, “Shaping the fluorescent emission by lattice resonances in plasmonic crystals of nanoantennas,” Phys. Rev. Lett. 102(14), 146807 (2009). [CrossRef]  

47. D. Punj, M. Mivelle, S. B. Moparthi, et al., “A plasmonic ‘antenna-in-box’ platform for enhanced single-molecule analysis at micromolar concentrations,” Nat. Nanotechnol. 8(7), 512–516 (2013). [CrossRef]  

48. K. L. Tsakmakidis, R. W. Boyd, E. Yablonovitch, et al., “Large spontaneous-emission enhancements in metallic nanostructures: towards LEDs faster than lasers [Invited],” Opt. Express 24(16), 17916–17927 (2016). [CrossRef]  

49. A. F. Koenderink, “Single-photon nanoantennas,” ACS Photonics 4(4), 710–722 (2017). [CrossRef]  

50. Z. Qi, Y. Zhong, and H. Liu, “Efficient method for the calculation of the optical force of multiple nanoparticles based on the coupling theory of quasinormal modes,” Opt. Lett. 46(18), 4610–4613 (2021). [CrossRef]  

51. H. Aouani, O. Mahboub, N. Bonod, et al., “Bright unidirectional fluorescence emission of molecules in a nanoaperture with plasmonic corrugations,” Nano Lett. 11(2), 637–644 (2011). [CrossRef]  

Supplementary Material (1)

NameDescription
Supplement 1       Supplement 1: includes Sections S1-S4 for providing some details of the derivation, the reciprocity-theorem method and the array scanning method, and includes Figs. S1-S4.

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (4)
Fig. 1.
Fig. 1. (a) Schematic diagram of the antenna of a single gold nanowire with the origin O located at the center of the antenna. (b) Distribution of Re(Ex) (top) and Re(Ey) (bottom) in z=0 plane for the electric-field x and y components of the QNM which is normalized to have Ex(L/2+10 nm,0,0)=1. (c1)-(c5) Distribution of |Ex| in z=0 plane at azimuthal angles ϕ=90°, 60°, 30°, 15°, and 0° for the scattered field of the antenna under the excitation of an x-polarized plane wave propagating along the negative z-direction. In the horizontal axis, r=|r| and the wavelength is λ=1127 nm. The red solid line, blue dash-dot line, and black dashed line represent the results obtained by full-wave FMM calculation, QNM expansion, and PESP, respectively. The insets show local enlargements within the range of r/λ<0.3.

View in Article | Download Full Size | PDF

Fig. 2.
Fig. 2. (a) Schematic diagram of a periodic array of single-nanowire antennas (array period a=b=1125 nm), where the red dot represents the point source (emission wavelength λ=1127 nm). (b) For an array size of 5×5, the far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$ calculated by the full-wave FMM (b1) and the QNM-coupling method (b2). (c) Positions of Rayleigh anomalies in the plane of $({\tilde{k}_x},{\tilde{k}_y})$, where the values of (m,n) in Eq. (13) are also provided. (d1)-(d3) $P({\tilde{k}_x},{\tilde{k}_y})$ calculated by the QNM-coupling method for array sizes of 15×15, 31×31 and 51×51, respectively. (e) For an infinite periodic array, $P({\tilde{k}_x},{\tilde{k}_y})$ calculated by the reciprocity-theorem method. (f1) Far-field emission-rate enhancement factor γrad=Γrad/Γvac calculated by the QNM-coupling method for different numbers of array periods (blue solid line), and the γrad for an infinite periodic array calculated by the reciprocity-theorem method (red dashed line). (f2) Total emission-rate enhancement factor γtotal=Γtotal/Γvac calculated by the QNM-coupling method for different numbers of array periods (blue solid line), and the γtotal for an infinite periodic array calculated by the ASM (red dashed line). The numbers of array periods along the x and y directions are set to be the same. (g) Distribution of the electric-field x-component Re(Ex) (normalized by Γvac) around the central antenna in the z=0 plane. The results are obtained by the QNM-coupling method for an array size of 51×51 (g1) and by the ASM for an infinite periodic array (g2).

View in Article | Download Full Size | PDF

Fig. 3.
Fig. 3. Far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$ for different values of the designated coupling distance RK. (a)-(h) RK=0 (i.e., setting QNM-coupling coefficients κ(p,m),(q,n)=0 for any p≠q), 5λ, 9λ, 15λ, 21λ, 31λ, 41λ, 51λ. The calculation is performed for a 51×51 periodic array of single-nanowire antennas, with parameters identical to those in Fig. 2(d3).

View in Article | Download Full Size | PDF

Fig. 4.
Fig. 4. (a) Schematic diagram of a periodic array of dipole nanoantennas, where the array period is a=b=1171 nm, and the point source (red dot) radiates at a wavelength λ=1171 nm. (b) Distributions in the z=0 plane for the electric field x-component Re(Ex) of symmetric QNM (upper image) and antisymmetric QNM (lower image) which are normalized to have Ex(d/2+L+10nm,0,0)=1. (c) Far-field radiation-intensity angular distribution $P({\tilde{k}_x},{\tilde{k}_y})$ for a 51×51 periodic array of antennas calculated with the QNM-coupling method. (d) $P({\tilde{k}_x},{\tilde{k}_y})$ for an infinite periodic array of antennas calculated with the reciprocity-theorem method using the full-wave FEM.

View in Article | Download Full Size | PDF

Tables (1)

Tables Icon

Table 1. Computation times of different methods. The computation is performed on a desktop PC with 6-core CPU of 4 GHz and a 64 GB RAM. The total number of sampling points of (θ,ϕ) for the reciprocity-theorem method in Fig. 2(e) and Fig. 4(d) is 1/4 of that for the QNM-coupling method. Therefore, the computation time given in the table for the reciprocity-theorem method is four times of the actual computation time.

View Table

Equations (15)

Equations on this page are rendered with MathJax. Learn more.

Ψ s ( r , ω ) = ∑ p = 1 P Ψ p s ( r , ω ) ,
Ψ p s ( r , ω ) = ∑ m = 1 M p α p , m ( ω ) Ψ ~ p , m ( r ) .
∇ × E ~ p , m = i ω ~ p , m μ 0 H ~ p , m , ∇ × H ~ p , m = − i ω ~ p , m ε p ( r , ω ~ p , m ) ⋅ E ~ p , m ,
[ ω I − K ( ω ) ] a ( ω ) = b ( ω ) .
κ ( p , m ) , ( q , n ) ( ω ) = − ω F p , m ∭ V p E ~ p , m ( r ) ⋅ Δ ε p ( r , ω ) ⋅ E ~ q , n ( r ) d 3 r
F p , m = ∭ R 3 { E ~ p , m ( r ) ⋅ ∂ [ ω ε p ( r , ω ) ] ∂ ω | ω = ω ~ p , m ⋅ E ~ p , m ( r ) − μ 0 H ~ p , m ( r ) ⋅ H ~ p , m ( r ) } d 3 r ,
β p , m ( ω ) = − ω F p , m ∭ V p E ~ p , m ( r ) ⋅ Δ ε p ( r , ω ) ⋅ E b ( r , ω ) d 3 r ,
Ψ p , μ s ( r , ω ) = ∑ m = 1 M p α p , m ( μ ) ( ω ) Ψ ~ p , m R ( r , ω ) ,
α p , m ( μ ) ( ω ) = [ − ω ∭ V p E ~ p , m ( r ) ⋅ Δ ε p ( r , ω ) ⋅ E p ( μ ) ( r , ω ) d 3 r ] / [ ( ω − ω ~ p , m ) F p , m ] ,
Ψ ~ p , m R ( r , ω ) = ∑ σ = 1 , 2 C ~ p , m + ( k ~ x , k ~ y , σ ) Ψ PW + ( k ~ x , k ~ y , σ ) exp ( i k 0 n b r ) λ i n b r ,
C ~ p , m + ( k ~ x , k ~ y , σ ) = k 0 2 n b η vac ∫ ∫ ∘ S ⁡ ( E PW − × H ~ p , m − E ~ p , m × H PW − ) ⋅ n ^ d s 8 π 2 [ E PW − ( − k ~ x , − k ~ y , σ ) ⋅ E PW + ( k ~ x , k ~ y , σ ) ] .
S ( k ~ x , k ~ y , r ) = ( λ n b r ) 2 ∑ σ = 1 , 2 | C ~ + ( k ~ x , k ~ y , σ ) | 2 S PW + ( k ~ x , k ~ y , σ ) ,
S PW + ( k ~ x , k ~ y , σ ) = 1 2 Re [ E PW + ( k ~ x , k ~ y , σ ) × H PW + ( k ~ x , k ~ y , σ ) ∗ ] ,
C ~ + ( k ~ x , k ~ y , σ ) = C ~ b ( k ~ x , k ~ y , σ ) + ∑ p = 1 P ∑ m = 1 M p α p , m C ~ p , m + ( k ~ x , k ~ y , σ ) g ( k ~ x , k ~ y , r p ) ,
( k x − m 2 π / a ) 2 + ( k y − n 2 π / b ) 2 = ( k 0 n b ) 2 ,
Select as filters
Select Topics Cancel
Top
       
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.
Privacy | Terms of Use