Resonantly enhanced second- and third-harmonic generation in dielectric nonlinear metasurfaces

Linear optical response of the metasurface

The efficiency of the frequency conversion processes in a nonlinear metasurface is chiefly determined by the distribution of the optical near-field around the meta-atoms and the strength of the light-matter interaction at the FF. Therefore, we start our investigations with the numerical calculation of linear transmittance, by sweeping the asymmetry parameter $\alpha$ and the angle of incidence $\theta$ and looking for resonances of the optical system as well as their linewidth and frequency dispersion properties. In these calculations, we used the rigorous coupled-wave analysis implemented in Synopsys’s DiffractMOD. The optical excitation is a plane-wave linearly polarized either along the long side of the bars ($x$-axis) or along their short side ($y$-axis).

In the case of $x$-polarized, normally incident plane-wave, as illustrated in Fig. 2(a), all spectral resonances are blue-shifted as the structural asymmetry parameter $\alpha$ increases from 0 to 1. When $\alpha$ becomes vanishingly small, the linewidth of the two sharp resonances at around 1780 nm vanishes so that their Q-factors diverge. This is a key feature of symmetry-protected BICs at the Γ point. Furthermore, based on the diffraction theory for 2D periodic structures, below the first Bragg-diffraction limit frequency, ${f}_{{\rm{B}}}=c/{n}_{{\rm{sub}}}\Lambda =240.8$ THz (${\lambda }_{{\rm{B}}}=c/{f}_{{\rm{B}}}=1245$ nm), where $c$ is the speed of light in vacuum and ${n}_{\rm{sub}}$ is the index of refraction of the substrate (glass), only the zeroth-order diffraction channel exists. For our design, the frequency of the two resonances is just below this diffraction limit, thus satisfying a key prerequisite for the existence of BICs. To utilize the physics of BICs, in-plane inversion symmetry is broken by reducing the width of one of the amorphous silicon bar pairs. The BICs generated via this procedure are expected to be characterized by a large Q-factor, chiefly due to the fact that radiative processes are suppressed.

Figure 2.  Simulated linear transmission spectra of the metasurface. (a) Dependence of transmission on the asymmetry parameter α under normal incidence and for x-polarized plane-wave excitation. The green dashed line corresponds to α = 0.075. (b) Dependence of transmission on the angle of incidence θ, determined for α = 0.075 and for x-polarized plane-wave excitation. (c, d) Same as in (a) and (b), respectively, but determined for y-polarized plane-wave excitation.

DownLoad: Full-Size Img PowerPoint

In addition to BICs, GMR-type resonances are also a promising way to realize metasurfaces that possess optical resonances with a large Q-factor. This is achieved through the coupling of diffracted waves propagating in the plane of the metasurface to the radiation continuum. For normal incidence, GMRs in two-dimensional periodic structures exist if the nanostructure satisfies the phase-matching condition for total internal reflection^{24, 56}: $\mathrm{max}\left[{n}_{{\rm{air}}},{n}_{{\rm{sub}}}\right]\le \lambda \left|\sqrt{{m}^{2}+{n}^{2}}/{Λ}\right| < {n}_{{\rm{eff}}}$, where ${n}_{{\rm{air}}}$ is the index of refraction of air, $\lambda$ is the working wavelength, $m$ and $n$ are integers defining the diffraction order, ${Λ}$ is the lattice period, and ${n}_{{\rm{eff}}}$ is the effective index of refraction of the unit cell. Using the value of ${n}_{{\rm{eff}}}$ calculated from the filling factor of the amorphous silicon in our designed metasurface⁵⁶, we determined that this critical condition defining the first diffraction order is satisfied in the wavelength range considered in this work, i.e., spanning the 1400 nm to 1900 nm domain. As a result, a GMR with a large Q-factor can be observed in the case of a symmetric meta-atom ($\alpha =0$) at 1622 nm under $x$-polarized normal incidence, as per Fig. 2(a). It should be noted that the linewidth of quasi-BICs is inversely proportional to the structural asymmetry parameter $\alpha$, whereas that of the GMR is practically independent on $\alpha$.

To gain a deeper understanding of the physical properties of BIC-type and GMR-type resonances, we further calculate the angle-dependent transmission map under the condition of $\alpha =0.075$, denoted as a green dashed line in Fig. 2(a), where the incident beam is polarized along the $x$-axis (s-polarized wave) and the angle of incidence is varied by changing the wavevector ${k}_{y}$ (${k}_{x}=0$). The corresponding results are shown in Fig. 2(b). Interestingly enough, the spectral locations of the two quasi-BICs diverge from each other as the angle of incidence increases, whereas the resonance wavelength of the GMR decreases monotonously as the angle of incidence $\theta$ increases up to about ${\theta \approx }$30°.

Due to the lack of ${C}_{4}$ rotational symmetry of the metasurface with bar-shaped elements in the unit cell, it is also relevant to study the linear optical response of the metasurface under s-polarized plane-wave excitation; the results of this analysis are summarized in Fig. 2(c) and 2(d). In the case when the meta-atom possesses in-plane inversion symmetry ($\alpha =0$), two BICs with vanishingly small linewidth are found at wavelengths of 1640 nm and 1736 nm, and a GMR with relatively large linewidth is found at 1760 nm. The evolution of these three resonances with the variation of $\theta$, at α = 0.075, is shown in Fig. 2(d), whereby the incident plane-wave is polarized along the $y$-axis ($s$-polarized plane-wave). The physical properties of the optical resonances under $y$-polarization of the incident plane-wave are similar to those in the $x$-polarized case, so that they will not be discussed in detail here. However, one fact that does need to be mentioned is that the resonant behavior of the quasi-BIC around 1715 nm at $\theta = 0^\circ$ in Fig. 2(d) weakens considerably when the angle of incidence increases beyond $2^\circ$. This finding can be explained by the large radiative losses into the free space when the incident wave is off Γ-point. We can thus conclude from the linear optical transmission maps for both polarizations that the metasurface possesses both BICs and high-Q GMRs. Importantly, the resonant behavior of these optical modes can be controlled by varying the asymmetry parameter and the angle of incidence of the plane-wave excitation.

The analysis of the modes of the metasurface can provide valuable information regarding the efficiency of the interaction between the incoming light and the metasurface. To this end, we considered the cubic lattice with asymmetric unit cell ($\alpha =0.075$) and calculated its eigenmodes using the finite-element method implemented in COMSOL Multiphysics. Due to the absence of in-plane inversion symmetry of the unit cell, here we investigate the band structure of the metasurface in the reciprocal lattice, along two paths in the $k$-space: $\left(0, 0\right)\to \left(0,{\rm{\pi}} /{Λ}\right)$ and $\left(0, 0\right)\to \left({\rm{\pi}} /{Λ}, 0\right)$. The results of this analysis are summarized in Fig. 3(a) and 3(b). In the band diagrams presented in these figures, only modes with Q-factors that are larger than 50 are presented, as only these modes would lead to a significantly enhanced light-matter interaction. In addition, this facilitates a comparison with the results revealed by the transmission spectral maps presented in Fig. 2.

Figure 3.  Eigenmode analysis in case of α = 0.075. Band diagrams with respect to (a) k_y and (b) k_x. Only modes with Q-factors larger than 50 are included. Red solid (blue dotted) lines are modes excited by x-polarized (y-polarized) plane-waves. (c–h) In-plane electric field profiles of the eigenmodes at the Γ point, presented in order of increasing wavelength. The yellow dashed lines indicate the interfaces between amorphous silicon and air.

DownLoad: Full-Size Img PowerPoint

The calculated eigenmodes present in the linear transmission spectra are classified according to the polarization of the incident plane-wave, which was either along $x$- or $y$-axis. As can be observed from Fig. 2(b) and Fig. 3(a), in the case of the eigenmodes that can couple to $x$-polarized plane waves, there is a good agreement between the frequency dispersion of resonances presented in the angle-resolved transmission map and in the eigenmode spectra (red solid lines). A similarly good agreement can be observed in the other case (eigenmodes that couple with $y$-polarized plane-waves), whereby the bands in the dispersion map in Fig. 2(d) agree well with the photonic bands indicated as blue dotted lines in Fig. 3(b). Moreover, Fig. 3(c–h) show the spatial distribution of the electric near-field ($x-y$ plane) within the unit cell in the middle cross-section of the amorphous silicon bar pairs at the Γ point. The field profiles correspond to the increasing wavelength of the optical eigenmodes. Due to the asymmetric nature of the unit cell of the dielectric metasurface, it is expected that the electric field patterns also present asymmetric features, which accounts for the optical coupling into free space of BIC-type resonances.

Since the rate at which a resonance of the optical metasurface radiates into the continuum is quantified by the Q-factor, it is particularly relevant to evaluate the Q-factor of the resonances. This physical quantity is important not just for the characterization of the light-matter interaction in the linear regime, but it also plays a critical role in the nonlinear optical regime whereby it is related to the optical near-field enhancement. We computed the dependence of the Q-factor at the Γ point of various resonances on the asymmetry parameter $\alpha$ of the meta-atom by using the eigenmode solver in COMSOL Multiphysics, the results of these calculations being shown in Fig. 4. It can be seen in this figure that for a symmetric metasurface ($\alpha =0$), for both polarizations of the incident plane-wave, there are eigenmodes with diverging Q-factors, implying that they are of a BICs nature. With the increase of the metasurface asymmetry parameter, the Q-factor of BICs decreases rapidly, which, as illustrated in Fig. 2(a) and 2(c), results in a considerable broadening of the corresponding linewidth. By contrast, the Q-factor of the GMRs is finite for all values of the asymmetry parameter and, as illustrated by the practically flat curves in Fig. 4, shows a weak dependence on the structural asymmetry parameter. Moreover, for the sake of completeness, for both BIC- and GMR-type resonances, we provide here the calculated resonance wavelengths and the corresponding Q-factors corresponding to the asymmetry parameter value ($\alpha =0.075$) at which the amorphous silicon metasurface was experimentally characterized: λ=1590 nm (Q = 2176, GMR), λ=1600 nm (Q = 400, BIC), λ=1715 nm (Q = 4723, BIC), λ=1747 nm (Q = 230, GMR), λ=1770 nm (Q = 1562, BIC), and λ=1778 nm (Q = 4965, BIC).

Figure 4.  Simulated dependence of Q-factor of the metasurface eigenmodes on the structural asymmetry parameter α for modes excited by (a) x-polarized and (b) y-polarized incident plane-waves. The calculations are performed at the Г point.

DownLoad: Full-Size Img PowerPoint

One key feature of BIC modes is that they have a topological nature, meaning that they are characterized by a topological invariant called topological charge. By using a well-known approach based on the dependence of the phase of the of the mode polarization in the $k$-space, we have found out that the topological charge of the BICs at 1600 nm, 1715 nm, 1770 nm, and 1778 nm is −1, 1, 1, and −1, respectively.

Enhanced SHG and THG from nonlinear metasurfaces

We now move on to the nonlinear optical response of the amorphous silicon metasurface and investigate the relationship between the excitation of BIC- or GMR-type resonances and the enhancement of nonlinear optical interactions. To illustrate the agreement between predictions of rigorous computational analysis and experimental results, we have fabricated an amorphous silicon metasurface corresponding to $\alpha =0.075$. Figure 5(a) presents a scanning electron microscope (SEM) image of the fabricated periodic pattern, showing the symmetry-broken feature of the square array dielectric metasurface. Before we characterized the nonlinear enhancement effects in the asymmetric metasurface, we measured the linear transmissivity spectra (see Fig. 5(b)) under $x$-polarized illumination, determined over a spectral interval ranging from 1500 nm to 2000 nm. As illustrated in Fig. 5(b), the transmission spectra obtained from experiments and simulation show similar features for wavelengths ranging from 1500 nm to 2000 nm. Particularly, a resonance with an asymmetric line shape can be observed at 1620 nm in the experimentally determined spectrum (marked by a purple arrow), corresponding to the sharp resonance seen at 1590 nm in the simulated transmission spectrum, the difference between the two values being less than 2%. This resonance is used to demonstrate further the dependence of the transmitted nonlinear optical signals at the SH and TH on the pump power at the FF. For the numerically calculated transmission spectrum, two more sharp resonances with extremely narrow linewidth exist at around 1780 nm. The fact that the signature of these extremely narrow resonances is not observed in the measured spectrum can be attributed to the inhomogeneous spectral broadening associated with inherent fabrication imperfections. Interestingly enough, as discussed above, these two resonances were found to be quasi-BICs.

Figure 5.  (a) Top-view SEM image of a fabricated metasurface corresponding to α = 0.075. (b) Comparison of experimentally measured linear transmittance and the corresponding numerically calculated spectrum at α = 0.075. The purple arrow denotes the resonance that is used to illustrate the measured power-dependent second-harmonic and third-harmonic signals presented in (c). (c) Measured power-dependent output SH (red dots) and TH (blue dots) signals from the fabricated sample with α = 0.075. The dependence of nonlinear emissions is measured at the resonance with (pump) wavelength of 1620 nm under field polarization oriented along the x-axis. The dots represent the measured data while the black dashed lines denote the fitting lines of SHG and THG power.

DownLoad: Full-Size Img PowerPoint

The dependence of the power generated at the SH and TH on the pump power, measured again for the metasurface with $\alpha =0.075$, is depicted in Fig. 5(c). As illustrated by this figure, the fitted slopes of the curves corresponding to SH and TH nonlinear optical interactions are 2.02 and 2.98, respectively, which clearly proves their quadratic and cubic dependence on the optical power at the FF, respectively (for pump power lower than 0.02 W, which was the power range where experiments were performed).

Considering the relation between the FF optical field and the SH/TH nonlinear polarizations described by Eqs. (1) and (2), one expects that the intensities of generated second- and third-harmonic are highly dependent on the enhancement of the local linear optical field. In addition, double-resonance mechanisms provide an effective method to further enhance nonlinear optical interactions. However, optical losses in amorphous silicon are large in the working spectral range of SHG (750 – 900 nm) and THG (500 – 600 nm), implying weak coupling between the resonances at the FF and the resonances in the vicinity of the SH/TH. Therefore, we do not provide here the linear transmission spectra at the SH and TH. In turn, we have calculated for both polarizations of plane-wave excitation at α = 0.075 the averaged enhancement of the linear local optical field, defined by the relation ${\eta }_{{\rm{E}}}=\dfrac{1}{V}{\displaystyle \int }_{V}\left|\boldsymbol{E}\left(\omega ;\boldsymbol{r}\right)\right|/{E}_{0}{\rm{d}}\boldsymbol{r}$, where ${E}_{0}$ is the amplitude of the incident electric field and $V$ is the volume of the bar-shaped resonators (amorphous silicon). The results regarding the enhancement of the local electric field are summarized in Fig. 6.

Figure 6.  Enhancement factor η_E of the averaged electric field amplitude within the amorphous silicon resonators under plane- wave illumination with the asymmetry of the metasurface being α = 0.075.

DownLoad: Full-Size Img PowerPoint

Comparing the field enhancement with the Q-factors in Fig. 4, we find that in case of $y$-polarization, resonances with higher Q-factor correspond to larger linear optical field enhancement. On the other hand, the resonances under $x$-polarized plane-wave excitation at 1590 nm and 1770 nm possess similar Q-factors and the corresponding enhancement of the linear field is of similar magnitude as well. Moreover, despite its high Q-factor, the enhancement corresponding to the resonance at 1778 nm ($x$-polarization) is relatively weak due to a smaller amount of field penetrating in the dielectric resonators of the metasurface (see field profile in Fig. 3(h)). Based on the above analysis, we are led to the conclusion that the second-order and third-order nonlinear optical interactions can be analyzed via the enhancement of the linear optical field, the Q-factor of the resonances, and the near-field distribution of the GMR and BIC resonances.

After establishing the nature of the nonlinear optical interactions, we turned our attention to a more detailed investigation of the fundamental physics pertaining to the generation of the SH signal. To this end, we first measured the spectral dependence of the SH intensity generated from the dielectric metasurface consisting of bar-shaped meta-atoms, in the case of $\alpha =0.075$. One key physical quantity that characterizes the spectrum of the generated SH signal is the SH enhancement factor, ${\eta }_{\mathrm{S}\mathrm{H}\mathrm{G}}$, defined as the ratio between the SH power emitted by the periodic array of dielectric resonators and the SH power generated from an amorphous silicon film of the same thickness. In the case of $x$-polarized plane-wave pump, three spectral peaks at wavelength of 1620 nm, 1684 nm, and 1714 nm at the fundamental wavelength can be observed in the SHG spectra [see the red curve in Fig. 7(a)]. Importantly, nearly 550 times enhancement of SHG intensity is achieved at the wavelength of 1620 nm, corresponding to the sharp resonance with a large Q-factor, marked by a purple arrow in Fig. 5(b). The additional two spectral resonances observed at longer wavelengths of 1684 and 1714 nm show smaller enhancement factors of 70 and 45, respectively.

Figure 7.  (a) Measured enhancement factor η_SHG of SHG from the all-dielectric metasurface with α = 0.075. An amorphous silicon slab with the same thickness as that of the metasurface is used as a reference. The inset shows a zoom-in of the enhancement factor. (b) Simulated linear transmittance for asymmetry parameter α = 0.075. (c, d) Numerically calculated surface and bulk contributions to SHG, respectively, when α = 0.075. Red (blue) lines correspond x-polarized (y-polarized)

DownLoad: Full-Size Img PowerPoint

To better understand the experimental results, we performed full-wave numerical simulations using the frequency solver of COMSOL Multiphysics and compared the computational results with the experimental data. To facilitate the identification of resonant features in the linear regime, we show in Fig. 7(b) the simulated transmissivity spectra for both polarizations at $\alpha = 0.075$, corresponding to the green dashed lines in Fig. 2(a) and 2(c), the calculations being performed for plane-waves normally incident onto the metasurface. For the nonlinear calculations, nonlinear optical interactions are described by the nonlinear surface and bulk polarizations given by Eq. (1). In this computational analysis, we separated the contributions to the SHG signal of the surface (local) and bulk (nonlocal) nonlinear polarizations, so that we unequivocally can establish which nonlinear optical effect dominates the nonlinear interactions. The results of this analysis are summarized in Fig. 7(c) and 7(d). It can be concluded from these figures that the main spectral features of the experimentally measured spectra are relatively well reproduced. For example, there is only a relative difference of about 2% between the measured and computed wavelengths of the main resonance located at 1620 nm in the measured SHG spectrum.

The comparison between experimental and simulation results reveals several additional ideas: First, the excitation of both the GMR (1620 nm) and quasi-BIC (1684 nm and 1714 nm) resonances produces a large enhancement of the generated SH. Second, the electric-dipole allowed nonlinear polarization is the dominant source of the SHG over almost the entire spectral domain considered in our study. This finding is consistent with the commonly used assumption that the bulk contribution to SHG can be neglected in the case of centrosymmetric materials. Third, there is a difference between the experimental and simulation values of the strength of SH enhancement corresponding to various resonances. These differences can be explained by the reduced experimental values of the Q-factors (as extracted from the experimental spectra of ${\eta }_{\mathrm{S}\mathrm{H}\mathrm{G}}$), a decrease attributable to structural perturbations such as fabrication imperfections and finite-size effects. Moreover, whereas in the experiments we measured the SH emitted in the zero-diffraction order, in numerical simulations we computed the total generated SH.

Importantly, despite the fact that the bulk SHG signal is negligible according to our numerically calculated results, the separate incorporation in the simulations of the two sources of SHG provides valuable information regarding the physical properties of the SHG interaction. For example, as per Fig. 7(c), the maximum SH intensity originating from surface effects is observed for the BIC-type resonance located at 1770 nm, whereas the largest value of SH intensity generated in the bulk is found to correspond to the GMR-type resonance located at 1590 nm. One can understand this finding by comparing the near-field distributions corresponding to the two resonances. With this in mind, one can conclude that it is possible to tune the ratio between the surface and bulk contributions to SHG from dielectric metasurfaces made of centrosymmetric optical materials, enabling new tools to use in the design of active photonic nanodevices.

We also measured SHG from the asymmetric metasurface under the condition of a $y$-polarized incident plane wave, the results being presented in Fig. 7(a). Comparing these results with the corresponding simulated spectra shown in Fig. 7(c) and 7(d), we see that again the spectral positions of the two ${\eta }_{\mathrm{S}\mathrm{H}\mathrm{G}}$ resonances (quasi-BIC at 1675 nm and GMR at 1700 nm) in the computed and measured spectra agree within a difference of about 2%. Overall, our simulated nonlinear spectra agree well with the measured SHG; however, perhaps due to the inhomogeneous spectral broadening effects, the resonance at 1600 nm in the computed spectra is not observed in experimentally measured ones. Using as a reference an amorphous silicon film with the same thickness as the metasurface, we determined that the SHG is resonantly enhanced by 110 times at 1675 nm and 70 times at 1700 nm. Similarly, to the case of $x$-polarized incident plane-waves discussed above, the contribution of the surface term to the SHG dominates when compared to the bulk contribution under $y$-polarized incident waves. Moreover, the large nonlinear enhancement arises chiefly from high-Q quasi-BIC at 1675 nm and the GMR at 1700 nm.

In addition to second-order nonlinear effects, the proposed resonant metasurface ($\alpha =0.075$) can also enhance THG arising from the large (bulk) third-order optical nonlinearity of amorphous silicon. To characterize the enhancement effect, we define a TH enhancement factor, ${\eta }_{\mathrm{T}\mathrm{H}\mathrm{G}}$, as the ratio between the TH power emitted by the asymmetric metasurface and the TH power generated from an amorphous silicon slab with the same thickness as that of the metasurface. The results of the measurements and the computed THG spectra corresponding to the third-order nonlinear polarization described by Eq. (2), determined under plane-wave excitation at the FH for both $x$-polarized (red) and $y$-polarized (blue) incident plane waves, are summarized in Fig. 8. Regarding the case of incident waves polarized along the $x$-axis, the high-Q GMR-type resonance at the FF wavelength (1620 nm) gives rise to extremely large THG enhancement, namely nearly 4200 times enhancement, whereas the two quasi-BICs originating from the symmetry breaking of the unit cell induce a ${\eta }_{\mathrm{T}\mathrm{H}\mathrm{G}}=340$ at 1674 nm and ${\eta }_{\mathrm{T}\mathrm{H}\mathrm{G}}=1341$ at 1715 nm. Moreover, the dielectric metasurface under an $y$-polarized incident wave excitation supports two resonances, namely a quasi-BIC at λ=1675 nm with enhancement factor of ${\eta }_{\mathrm{T}\mathrm{H}\mathrm{G}}=3790$ and a GMR at λ=1712 nm with ${\eta }_{\mathrm{T}\mathrm{H}\mathrm{G}}=4930$.

Comparing Fig. 8(a) and 8(b) we see that, overall, for both polarizations there is a good agreement between the experimental results and numerically computed spectra characterizing the third-order nonlinear optical response of the asymmetric bar-shaped metasurface. Due to fabrication and structural imperfections, there is a small difference between the resonance wavelengths of the peaks seen in the measured and computed nonlinear enhancement spectra. In addition, due to inhomogeneous spectral broadening, the Q-factor of measured spectral resonances is somewhat smaller than that of the computed spectra (see the quasi-BIC at λ=1675 nm excited by an $y$-polarized plane-wave). This effect results in a decrease of THG intensity. It should also be noted that despite the fact that the computed spectra for $x$-polarized incident plane waves display two quasi-BIC resonances at $\lambda \approx 1770\mathrm{ }\mathrm{n}\mathrm{m}$, only one such resonance can be observed experimentally at this frequency. This is so because the two corresponding resonances at the FF merge as the result of inhomogeneous spectral broadening effects induced by fabrication imperfections.

Figure 8.  (a) Measured enhancement factor of THG from the resonant metasurface described by α = 0.075. (b) Numerically calculated TH intensity spectra corresponding to the same metasurface. Red (blue) lines correspond to x-polarized (y- polarized) incident plane-waves. In both figures, the upper ticks of the x-axis denote the TH wavelength.

DownLoad: Full-Size Img PowerPoint

According to our experimental measurements, the up-conversion efficiencies of second-order (${P}_{\mathrm{S}\mathrm{H}\mathrm{G}}/{P}_{\mathrm{p}\mathrm{u}\mathrm{m}\mathrm{p}}$) and third-order (${P}_{\mathrm{T}\mathrm{H}\mathrm{G}}/{P}_{\mathrm{p}\mathrm{u}\mathrm{m}\mathrm{p}}$) nonlinear processes are estimated to be ${10}^{-8}$ and ${10}^{-7}$, respectively. To compare these values to the state-of-the-art, we mention that the current record efficiencies are $6.8\times {10}^{-7}$ for SHG from arrayed silicon slotted nanocubes⁴¹ and $1.8\times {10}^{-6}$ for THG measured from silicon pillar metasurfaces⁵⁷. Thus, the nonlinear conversion efficiency reported in this work is not record high but comparable to the state-of-the-art.

It is also relevant to consider the directionality of the SHG and THG transmission, since the working wavelengths at SHG and THG are beyond the diffraction limit (1245 nm). To characterize the nonlinear emissions into different diffraction channels, we performed the Fourier analysis of numerically calculated transmitted nonlinear signals, for the cases of SHG at 1620 nm under $x$-polarized excitation, which corresponds to the highest measured SHG enhancement, and THG at 1747 nm under $y$-polarized excitation, which corresponds to the largest measured THG enhancement.

Figure 9 displays the normalized power distribution in each diffraction channel, where $m$ and $n$ denote the orders of diffraction corresponding to the $x$- and $y$-directions, respectively. Regarding the SHG (see Fig. 9(a)), our results indicate that 9 diffraction channels exist and the emission angle ${\theta }_{\mathrm{S}\mathrm{H}\mathrm{G}}$ for each of them are: 0° for (0,0) order; 39.7° for (1, 0), (−1, 0), (0, 1) and (0, −1) orders; and 64.6° for (1, 1), (1, −1), (−1, 1) and (−1, -1) orders. Most of the nonlinear power at SH radiates into the (1, 0) and (−1, 0) channels, whereas the emission into zeroth order is only about $0.05$ of the power emitted in the diffraction order corresponding to the highest emission. Regarding THG, the periodic structure of the metasurface allows 13 diffraction orders at the TH wavelength, and the emission angle ${\theta }_{\mathrm{T}\mathrm{H}\mathrm{G}}$ for these diffraction channels are: 0° for (0, 0) order; 27.9° for (1, 0), (−1, 0), (0, 1) and (0, −1) orders; 41.4° for (1, 1), (1, −1), (−1, 1) and (−1, −1) orders; and 69.3° for (2, 0), (−2, 0), (0, 2) and (0, −2) orders. Clearly, most power is emitted in the (0, −2) diffraction order, and compared to this, the THG emission in the zeroth order is about $0.14$. All diffraction orders are included in our simulations by integrating the outgoing power flow in the transmission side of the metasurface. In the case of experimental measurements, the nonlinear signals are collected only from the zero-diffraction order.

Importantly from a practical point of view, we have demonstrated in this section both theoretically and experimentally that strong enhancement of several nonlinear optical processes (SHG and THG in our case) can be achieved in the same nonlinear metasurface, which suggests that such photonic nanodevices can become instrumental in designing active optical devices with new or improved functionality.

Dependence of SHG and THG on the metasurface asymmetry

As discussed above from both experimental and computational perspectives, high-Q resonances supported by dielectric metasurfaces make it possible to realize enhanced second- and third-order nonlinear optical interactions. The calculated transmission spectra exhibited in Fig. 2 demonstrate that the linewidth (or Q-factor) of both GMR-type and BIC-type resonances, which determines the local field enhancement, are sensitive to the variation of the asymmetry of the meta-atom. We now go a step further aiming to gain a deeper and more comprehensive insight into the dependence of harmonic generation on the asymmetry features of the dielectric metasurface. To this end, we simulate the nonlinear optical response of the metasurface for different values of $\alpha$, that is for 0, 0.075, 0.15, and 0.225. Since we observed a good agreement between experimental and simulation data, here we present only results obtained through rigorous simulations of the bar-shaped amorphous silicon metasurface.

Figure 9.  Fourier analysis of numerically calculated (a) transmitted SHG at 1590 nm (under x-polarized incident plane-wave), which corresponds to the resonant peak of highest measured SHG enhancement, and (b) transmitted THG at 1747 nm (under y- polarized incident plane-wave), which corresponds to the resonant peak of largest THG enhancement in experiments. Fourier analysis is performed for an asymmetry parameter of α = 0.075.

DownLoad: Full-Size Img PowerPoint

The simulated SHG spectra presented in Fig. 10 reveal that the resonance wavelength and the linewidth of the corresponding resonances vary considerably with the asymmetry parameter, $\alpha$. Comparing the two contributions to SHG, we can infer that for both polarizations and for all values of $\alpha$, the surface contribution to the SHG is much larger than that of the bulk term. It is worth noting that in the case of $\alpha =0$, BICs are completely uncoupled with the radiation continuum and manifest themselves as spectral peaks with vanishingly small width [see the red lines depicting the SHG intensity spectra in Fig. 10(a–d)]. As $\alpha$ increases, due to the introduced structural asymmetry, the BICs are no longer symmetry-protected, which results in a strong optical coupling between the localized modes and the extended propagating waves. The signature of this strong optical interaction is the appearance of sharp resonances in the nonlinear response spectra. Furthermore, the results of the eigenmode analysis presented in Fig. 4 show that the Q-factors of BIC-type resonances decrease with the asymmetry parameter. Importantly, an opposite trend is observed in the case of the GMR-type resonances, namely their Q-factors slowly increase when the asymmetry parameter increases. It is expected therefore that the corresponding surface and bulk SHG present a similar dependence on the structural asymmetry. To confirm this conjecture, for both polarizations we analyze the dependence of each resonance of the metasurface on the asymmetry parameter $\alpha$. As can be seen in Fig. 10(a–d), indeed the width of the SHG resonances corresponding to quasi-BICs gradually increases with the increase of asymmetry parameter, whereas the opposite is true in the case of GMR resonances.

Finally, we investigated the relationship between the asymmetry parameter and the efficiency of the THG process. The corresponding simulated results are plotted in Fig. 11, where a strong dependence of TH intensity on the asymmetry parameter $\alpha$ can be clearly observed for plane-wave excitation at the FF wavelength with the field polarized along either $x$- or $y$-axis. Comparing these results with the SHG spectra in Fig. 10, one can draw conclusions similar to those valid in the case of SHG: First, both quasi-BICs and GMR-type resonances with large Q-factors induce greatly enhanced TH nonlinear signals. Second, the spectral location of resonances at the TH intensity, corresponding to different values of the asymmetry parameter $\alpha$, are closely related to the resonances at the FF. Third, the width of BIC-induced resonances at TH increases with the asymmetry of the metasurface, whereas the opposite trend is observed at the resonance wavelengths corresponding to excitation of GMR-type resonances. We also note that in the experimental samples both SHG and THG efficiency dependence on the asymmetry is affected by the ratio of the radiative and other types of losses. The maximal enhancement can be usually achieved in the critical coupling regime⁴⁰.

Figure 10.  Simulated dependence of SHG from the amorphous silicon metasurface with respect to the structural asymmetry parameters α in cases of (a) surface contribution under x-polaried excitation, (b) surface contribution under y-polarized excitation, (c) bulk contribution under x-polarized excitation, and (d) bulk contribution under y-polarized excitation.

DownLoad: Full-Size Img PowerPoint

Figure 11.  Numerical simulation of asymmetry-dependent TH intensity under plane wave excitation with polarization along (a) x and (b) y directions, respectively. In both figures, the upper ticks of the x-axis denote the TH wavelength.

DownLoad: Full-Size Img PowerPoint