1. Introduction

Metasurfaces, which are artificial optical thin films structured at subwavelength scale, have become widely used tools to engineer in unique ways light-matter interaction at the nanoscale [1–5]. These developments have been primarily facilitated by recent advancements in nanofabrication techniques, and have made possible to achieve planar structures that resonantly interact with light at specific frequencies. Compared to three-dimensional (3D) metamaterials, two-dimensional (2D) metasurfaces allow the manipulation of light using more compact devices [6] and exhibit lower optical power loss [7]. In addition, key nonlinear optical functionalities can be more easily implemented using metasurfaces [8–11], since usually one does not have to be concerned with satisfying phase-matching conditions. These and other unique optical properties of metasurfaces have been opening up new avenues to pursue promising applications in areas such as quantum information processing [12,13] and photonic integrated circuits[14, 15].

There are two main approaches towards nonlinear metasurfaces, namely metasurfaces based on metallic (plasmonic) [8, 10, 16–18] and dielectric [4, 19–21] materials. Plasmonic metasurfaces can generate strong optical near-field enhancement, a key feature for nonlinear optics applications. This comes, however, at the price of relatively large optical losses, which reduces their conversion efficiency. Dielectric metasurfaces, on the other hand, are characterized by small optical losses but their limitations stem from the fact that the enhancement of the optical near-field is smaller than that achieved with plasmonic metasurfaces. In order to overcome this shortcoming of all-dielectric metasurfaces, one could engineer their primary building blocks so that they possess resonant optical modes with high quality-(Q) factor. These specially engineered modes are different in nature from Mieresonances of dielectric particles [22, 23], which are characterized by relatively small Q-factor. Due to the large set of parameters that can be tuned to optimize their Q-factor, photonic crystal (PhC) cavities represent an excellent choice for designing optical resonators with high Q-factor, although other solutions, such as whispering-gallery modes of resonators made of silica or other dielectric materials, are possible, too.

In order to illustrate this idea, in this paper we demonstrate that the effective Raman susceptibility of a metasurface consisting of a periodic array of PhC cavities in a silicon PhC slab waveguide can be enhanced by more than 4 orders of magnitude, as compared to the intrinsic Raman susceptibility of silicon. This is achieved by designing a silicon PhC cavity that possesses two high-Q optical modes separated spectrally by the Raman frequency of silicon, Ω/2π = 15.6 THz [24]. The effective Raman susceptibility of the metasurface is calculated using a new homogenization procedure based on the effective-medium theory, an approach that is valid beyond the subwavelength unit cell regime. Importantly, this procedure can be applied to other multi-wavelength nonlinear optical interactions, such as four-wave mixing, optical parametric amplification, and second-harmonic generation, simply by designing PhC cavities that possess optical modes at the required frequencies. The paper is organized as follows. In the next section we present the geometrical configuration and optical properties of the nonlinear metasurface. Then, in Sec. 3 we introduce the homogenization method used to calculate the effective Raman susceptibility of the metasurface, whereas in Sec. 4 we present and discuss the main results. Finally, in Sec. 5 the main conclusions are summarizad.

2. Geometrical configuration and properties of the nonlinear metasurface

The nonlinear photonic metasurface and the homogenization procedure used to compute the corresponding effective Raman susceptibility are schematically illustrated in Fig. 1. It comprises of a rectangular array of PhC cavities in a PhC slab waveguide made of silicon (n_Si = 3.4). The slab waveguide consists of a 2D hexagonal lattice of air holes in a silicon slab, with lattice constant, a, hole radius, r = 0.29a, and slab thickness, t = 0.6a. Moreover, in this study we use L5 PhC cavities, namely they are formed by filling in 5 consecutive holes located on a line oriented along the ΓK symmetry axis of the lattice. The center-to-center distance between adjacent PhC cavities along the longitudinal (x) and transverse (y) directions is d_l = 17a and $d_t=6\sqrt{3}a$, respectively. In order to increase the Q-factor of the optical modes of the cavity the end holes are shifted by s_l = 0.15a [25].

 

Fig. 1 Schematic of a silicon nonlinear metasurface for enhancement of stimulated Raman scattering. The metasurface consists of a rectangular array of PhC cavities in a hexagonal PhC slab waveguide made of silicon. Also illustrated is the homogenization procedure used to compute the corresponding effective Raman susceptibility, ${\mathbf{\chi }}_R^{(3),\text{eff}}$, of a homogeneous slab with the same thickness, t, from the intrinsic Raman susceptibility of silicon, ${\mathbf{\chi }}_R^{(3)}$.

Download Full Size | PDF

The photonic band structure of the PhC slab waveguide and the frequencies of the cavity modes were computed using RSoft’s BandSOLVE [26] and are presented in Fig. 2, whereas the Q-factors of the cavity modes were determined using MEEP [27], a freely available software that implements the finite-difference time-domain algorithm. Thus, the PhC cavity possesses two optical modes with normalized frequencies (measured in units of 2πc/a) of ${\overline{\omega }}_p=0.2778$ and ${\overline{\omega }}_p=0.2605$ and Q-factors of Q_p = 1804 and Q_S = 1.12 × 10⁵. These cavity modes lie in the transverse-magnetic band-gap of the PhC slab waveguide, as per Fig. 2(b), the field profiles of the two optical modes being shown in this figure, too. Moreover, if the lattice constant is a = 333 nm, the frequencies of the two cavity modes are ω_p = 1572.5 THz and ω_S = 1474.6 THz, that is they are separated by the Raman frequency of silicon, ω_p − ω_S = Ω = 2π × 15.6 THz. Under these circumstances, the PhC cavity can be viewed as a “meta-molecule” possessing strong Raman response, as schematically illustrated in Fig. 2(a).

 

Fig. 2 (a) Schematic diagram of stimulated Raman scattering: a pump photon with frequency ω_p interacts with a Raman-active optical medium and generates a Stokes photon with frequency ω_S and a phonon with frequency Ω = ω_p − ω_S. (b) Transverse-magnetic band structure of the photonic crystal and the two cavity modes. The field profile of the cavity modes are shown in the right panels.

Download Full Size | PDF

Consider now that such PhC cavities are arranged in a rectangular 2D array, such that one creates a nonlinear metasurface that can potentially generate strongly enhanced stimulated Raman scattering as compared to that of intrinsic silicon. This is so because the existence of the two cavity modes ensures an effective optical coupling between the incident optical wave and the cavity mode at the pump frequency and provides an efficient mechanism for the Raman signal to be radiated in the continuum via the excitation of the cavity mode at the Stokes frequency. We have computed the effective Raman susceptibility of such a nonlinear optical metasurface using a newly developed homogenization technique, which we present in the next section.

3. Homogenization approach for calculation of the effective Raman susceptibility of the nonlinear metasurface

The most common approach to homogenize a metasurface is based on the effective-medium theory [28–31]. Generally speaking, this approach amounts to averaging the electromagnetic field over the unit cells of the metasurface. Through the homogenization procedure, the photonic system is reduced to a slab of homogeneous optical medium characterized by effective electric permittivity and magnetic permeability [32–34]. On the other hand, in the nonlinear case, the nonlinear effective susceptibilities are calculated form a spatial overlap integral among fields excited in the structure by probing waves at pump and higher-harmonics, weighted by the local nonlinear susceptibility tensor and averaged over the volume of the unit cell [35–38].

Before describing our homogenization method, we would like to point out that whereas the linear effective optical constants, such as the electric permittivity and magnetic permeability, can be calculated in an unambiguous way, in the case of the effective nonlinear susceptibilities the large number of nonvanishing components of these tensor quantities could lead to certain ambiguity in how they are defined and calculated. To be more specific, in most cases homogenization methods for nonlinear metasurfaces rely on the condition that the averaged nonlinear polarization of the metasurface is equal to the nonlinear polarization in the nonlinear homogenized slab of material. As mentioned, this approach results in ambiguities in the determination of the effective nonlinear susceptibility because three relations that ensure that the two nonlinear polarizations are identical must be used to determine a much larger number of components of the nonlinear susceptibility tensor. The key novel aspect of our homogenization approach is that the calculated effective nonlinear susceptibility is unique. This is achieved by imposing the condition that the averaged nonlinear polarization of the metasurface and the nonlinear polarization in the homogenized slab are term-wise identical when expressed in terms of the components of the nonlinear susceptibilities and the electric field components.

The Raman nonlinear optical response of the metasurface is determined by the nonlinear Raman polarization at the Stokes frequency, ω_S [39]:

The position-dependent Raman polarization defined by Eq. (1) can be used to calculate the spatially averaged effective Raman polarization:

In this equation and in what follows, the bar symbol indicates that the field corresponds to the homogenized medium. Moreover, we stress that in Eq. (3) the fields ${\overline{\mathbf{E}}}_{\text{eff}}({\omega }_p)$ and ${\overline{\mathbf{E}}}_{\text{eff}}({\omega }_S)$ are the effective optical fields at the pump and Stokes frequencies, respectively, generated by an input field that is equal to the field that excites the metasurface. They are given by:

The nonlinear optical response of the metasurface and homogenized slab are identical if the effective Raman polarizations described by Eqs. (2) and (3) are equal, which leads to three equations (one for each component). Note, however, that in the general case the effective Raman susceptibility tensor, ${\mathbf{\chi }}_R^{(3),\text{eff}}$, has 81 independent components, so that the corresponding system of equations is overdetermined. In order to circumvent this issue, we impose the condition that the r.h.s. of Eqs. (2) and (3) are term-wise identical. Using this constraint, it can be seen that the components of the effective Raman susceptibility tensor are determined by the following relations:

Note that the components of ${\mathbf{\chi }}_R^{(3),\text{eff}}$ and ${\mathbf{\chi }}_R^{(3)}$ cancel for the same set of indices i, j, k, and l.

The remaining physical quantities needed to calculate the effective Raman susceptibility given in Eq. (5) are the effective fields in the homogenized slab. They can be calculated analytically (see the Appendix), but before this one needs to determine the effective electric permittivity of the slab. This can be calculated by imposing the condition that the spatially averaged electric displacement in the metasurface, $D_i({\omega }_{\alpha })=(1/V){\displaystyle {\sum }_j{\displaystyle {\int }_V{ϵ}_{ij}(\mathbf{r})E_j({\omega }_{\alpha };\mathbf{r})d\mathbf{r}}}$, where ${ϵ}_{ij}(\mathbf{r})={ϵ}_0{n}_{Si}^2{\delta }_{ij}$ and ϵ_ij (r) = ϵ₀δ_ij for r in the silicon and air regions, respectively, with δ_ij being the Kronecker delta, is equal to the effective electric displacement in the homogenized slab, ${\overline{D}}_i({\omega }_{\alpha })=(1/V){\displaystyle {\sum }_j{\overline{ϵ}}_{ij}^{\text{eff}}({\omega }_{\alpha }){\displaystyle {\int }_V{E}_j({\omega }_{\alpha };\mathbf{r})d\mathbf{r}}}$. If the two electric displacement fields are equal on a term by term basis, the effective permittivity is given by:

This equation implies that the effective electric permittivity of the homogenized slab is described by a diagonal matrix, but the diagonal elements are not necessarily the same.

Finally, the effective permittivity in Eq. (6) is used to determine analytically the effective electric field in the homogenized slab. The details of these calculations are provided in the Appendix, so that here we only give the final result:

Here, E_i,0, i = x, y, is the amplitude of the incident field and

Before we move on to illustrate how this homogenization method can be applied to our metasurface, we stress that although the method can be successfully used to describe metasufaces with period comparable or larger than the wavelength of the interacting beams, the resulting effective permittivity and effective nonlinear Raman susceptibility would depend on the polarization of the interacting beams and the angle of incidence. Thus, it is known that when the wavelength is comparable to the characteristic length of a photonic structure, e.g. the lattice constant of a PhC, physical quantities related to the spatial average of the microscopic electromagnetic field depend not only on the frequency but also on the wave vector, k. This means that, for example, the dispersion relation of the electric permittivity has a more complex form, ϵ = ϵ(ω, k). This is valid for the nonlinear optical response of photonic structures, too, so that the nonlinear susceptibilities depend on the direction of propagation and polarization of the interacting fields. This is expected because when the wavelength is comparable to the characteristic length of a photonic structure the local fields, which determine the nonlinear polarization of the medium, depend on the direction of propagation and polarization of the incident fields.

4. Results and discussion

We have applied the homogenization procedure described in the preceding section to our nonlinear metasurface. The first step in our analysis was to compute the field distribution at the pump and signal frequencies of the corresponding cavity modes under plane-wave excitation conditions. The calculations were performed using the frequency-domain finite element method implemented in CST Studio, a commercially available software [41]. Thus, we first determined the spectral response of the metasurface and then, using a normally incident plane-wave excitation source polarized at π/4 with respect to the x-axis and with excitation frequency equal to the pump and Stokes frequencies we determined the field profiles at the corresponding resonance frequencies. The spectral response of the metasurface has been determined in the frequency domain ranging from 200 THz to 300 THz, that is a spectral domain that contains the two resonances of interest, by determining the frequency dependence of the electric field at an arbitrary point in the metasurface. Despite the fact that the spectrum varies with the position of the point in the metasurface where the electric field is probed, the location of the resonances does not depend on this choice.

We present in Fig. 3 the spectrum calculated using this approach, as well as the profiles of the field amplitude corresponding to the two cavity modes located at the pump and Stokes frequencies. These are the same modes shown in the inset of Fig. 2(b) and, as mentioned, were calculated using BandSOLVE. It can be seen that the spectral separation between the two resonances is equal to the Raman frequency, Δν = ν_R = 15.6 THz, a result that validates our design approach.

 

Fig. 3 Spectrum of the optical response of the metasurface.

Download Full Size | PDF

Due to the symmetry properties of the metasurface and the orientation of the cavity array with respect to the principal axes of silicon, the only non-zero component of the effective Raman susceptibility is ${\chi }_{R,1212}^{(3),\text{eff}}$. Therefore, in order to compute it, one only needs to calculate the effective permittivity components ${\overline{ϵ}}_{ij}^{\text{eff}}$, i, j = x, y. Note that the effective permittivity must be calculated both at the pump and Stokes frequencies. The effective permittivity is determined using Eq. (6), the results being ${\overline{ϵ}}_{xx}({\omega }_s)=8.777+0.367\mathrm{i}$, ${\overline{ϵ}}_{yy}({\omega }_s)=9.034-0.003\mathrm{i}$, ${\overline{ϵ}}_{xx}({\omega }_p)=9.000+0.044\mathrm{i}$, ${\overline{ϵ}}_{yy}({\omega }_p)=8.921+0.049\mathrm{i}$, and ${\overline{ϵ}}_{xy}({\omega }_s)={\overline{ϵ}}_{yx}({\omega }_s)={\overline{ϵ}}_{xy}({\omega }_p)={\overline{ϵ}}_{yx}({\omega }_p)=0$. As expected, the effective permittivity is a complex, diagonal matrix, and the diagonal elements are not identical.

The last step of the homogenization procedure is to use Eq. (5) to calculate the effective Raman susceptibility of the metasurface, ${\chi }_{R,1212}^{(3),\text{eff}}$. For the nonlinear metasurface made of L5 PhC cavities described in Fig. 1 our calculations show that ${\chi }_{R,1212}^{(3),\text{eff}}=-i18.8\times {10}^{-14}{\mathrm{m}}^2{\mathrm{V}}^{-2}$

, that is the effective Raman susceptibility is enhanced by $\eta =|{\chi }_{R,1212}^{(3),\text{eff}}/{\chi }_{R,1212}^{(3)}|=1.67\times {10}^4$ as compared to that of silicon. Interestingly enough, similar values of the enhancement of the Raman interaction in silicon photonic structures have been observed in the case of silicon PhC waveguides operated in the slow-light regime [42].

In order to understand the dependence of the effective Raman susceptibility on the size of the unit cell of the metasurface, we varied the number of holes, N, located along the x-direction between the end of the cavity and the boundary of the unit cell. As illustrated in Fig. 4, in which we plot this dependence, the effective Raman susceptibility of the nonlinear metasurface decreases when N increases and remains practically unchanged for N ≥ 4.5. The maximum enhancement is η = 4.22 × 10⁴, that is an additional increase by a factor of about 3. This dependence is readily explained by the fact that as N decreases there are more cavities per unit area and consequently the nonlinear optical response of the metasurface increases. It should be noted, however, that as the cavities are more closely packed together the mutual optical coupling can lead to a frequency shift of the resonance frequencies of the modes and thus to a decrease of the efficiency of the Raman interaction. This can be avoided by simply designing the cavity such that the frequency difference between the modes differs from the Raman frequency by exactly the frequency shift induced by the optical coupling between the cavities.

 

Fig. 4 The enhancement of the effective Raman susceptibility of the nonlinear metasurface relative to the intrinsic Raman susceptibility of silicon vs. the number of holes located along the x-direction between the end of the nanocavity and the boundary of the unit cell.

Download Full Size | PDF

5. Conclusions

In conclusion, by designing a nonlinear optical metasurface made of a two-dimensional periodic array of silicon photonic crystal cavities, we have demonstrated that the effective nonlinear Raman susceptibility of the metasurface can be enhanced by more than 4 orders of magnitude as compared to that of silicon. In order to achieve this dramatic enhancement of the nonlinear optical response of the metasurface, the photonic crystal cavities are designed so as to possess two optical cavity modes spectrally separated by the Raman frequency of silicon. Importantly, the linear and nonlinear optical response of the metasurface have been quantified using a new homogenization method, which is valid even when the period of the metasurface is comparable or larger than the wavelength of the interacting beams. The ideas presented in this paper have wide applicability, as they can be easily extended to other nonlinear optical interactions of practical interest, including second-and third-harmonic generation, four-wave mixing, and sum-and difference-frequency generation.

Appendix: Average electric field in a homogeneous slab of material

In this Appendix we present a derivation for the averaged electric field in a slab of homogeneous material. The corresponding electromagnetic structure and field configuration are presented in Fig. 5. Thus, we assume that a homogeneous slab, Medium 2, with permittivity ε₂ is sandwiched in-between two media, Medium 1 and Medium 3, which have permittivities ε₁ and ε₃, respectively. For simplicity, we assume that all media have the same permeability, µ₁ = µ₂ = µ₃ = µ₀. Moreover, incident on the slab is a plane wave propagating along the z-axis, with the electric and magnetic fields written as ${\mathbf{E}}_{1i}(z)=\widehat{\mathbf{i}}{E}_{1i}{e}^{i{k}_{1}z}$ and ${\mathbf{H}}_{1i}(z)=\widehat{\mathbf{j}}(E_{1i}/Z_1)e^{i{k}_{1}z}$, respectively, where $Z_1=\sqrt{{\mu }_0/{ϵ}_1}$ and $k_1=\omega \sqrt{{\mu }_0{ϵ}_1}$. Note that in this derivation we assume that ε₂ is diagonal but the diagonal elements are not necessarily equal. Under these circumstances, the fields in the three regions can be written as:

 

Fig. 5 Wave configuration corresponding to a homogeneous slab (Medium 2) with thickness, t, sandwiched in-between two media (Medium 1 and Medium 3).

Download Full Size | PDF

Medium 1

(9)$${\mathbf{E}}_1(z)={\mathbf{E}}_{1i}(z)+{\mathbf{E}}_{1r}(z)=\widehat{\mathbf{i}}\left(E_{1i}{e}^{i{k}_{1}z}+E_{1r}{e}^{-i{k}_{1}z}\right),$$

(10)$${\mathbf{H}}_1(z)={\mathbf{H}}_{1i}(z)+{\mathbf{H}}_{1r}(z)=\widehat{\mathbf{j}}\frac{1}{Z_1}\left(E_{1i}{e}^{i{k}_{1}z}-E_{1r}{e}^{-i{k}_{1}z}\right).$$

Medium 2

(11)$${\mathbf{E}}_2(z)={\mathbf{E}}_{2+}(z)+{\mathbf{E}}_{2-}(z)=\widehat{\mathbf{i}}\left(E_{2+}{e}^{i{k}_2(z-t}+E_{2-}{e}^{-i{k}_2(z-t)}\right),$$

(12)$${\mathbf{H}}_2(z)={\mathbf{H}}_{1+}(z)+{\mathbf{H}}_{1-}(z)=\widehat{\mathbf{j}}\frac{1}{Z_2}\left(E_{2+}{e}^{i{k}_2(z-t}-E_{2-}{e}^{-i{k}_2(z-t)}\right).$$

Medium 3

(13)$${\mathbf{E}}_3(z)={\mathbf{E}}_{3t}(z)=\widehat{\mathbf{i}}{E}_{3t}{e}^{i{k}_{3}z},$$

(14)$${\mathbf{H}}_3(z)={\mathbf{H}}_{3t}(z)=\widehat{\mathbf{j}}\frac{E_{3t}}{Z_3}{e}^{i{k}_{3}z}.$$

To simplify the notations, let us introduce the reflection coefficients at z = 0, Γ(0) = E_1r/E_1i, and at z = t, Γ(t) = E₂₋/E₂₊. Then, the fields in Medium 1 and Medium 2 can be expressed as:

(15)$${\mathbf{E}}_1(z)=\widehat{\mathbf{i}}{E}_{1i}\left(e^{i{k}_{1}z}+\mathrm{\Gamma }(0)e^{-i{k}_{1}z}\right),$$

(16)$${\mathbf{H}}_1(z)=\widehat{\mathbf{j}}\frac{E_{1i}}{Z_1}\left(e^{i{k}_{1}z}-\mathrm{\Gamma }(0)e^{-i{k}_{1}z}\right),$$

(17)$${\mathbf{E}}_2(z)=\widehat{\mathbf{i}}{E}_{2+}\left(e^{i{k}_2(z-t)}+\mathrm{\Gamma }(t)e^{-i{k}_2(z-t)}\right),$$

(18)$${\mathbf{H}}_2(z)=\widehat{\mathbf{j}}\frac{E_{2+}}{Z_2}\left(e^{i{k}_2(z-t)}-\mathrm{\Gamma }(t)e^{-i{k}_2(z-t)}\right).$$

The continuity of the tangent components of the fields at z = 0 and z = t requires that:

(19)$$\mathrm{\Gamma }(0)=\frac{Z_2(Z_3-Z_1)-i(Z_2^2-Z_1{Z}_3)\tan (k_{2}t)}{Z_2(Z_3+Z_1)-i(Z_2^2+Z_1{Z}_3)\tan (k_{2}t)};\mathrm{\Gamma }(t)=\frac{Z_3-Z_2}{Z_3+Z_2}.$$

In addition, using Eqs. (9)–(12) and the continuity of the electric field at z = 0, one can fiend the amplitude of the electric field in Medium 2

(20)$$E_{2+}=\frac{1+\mathrm{\Gamma }(0)}{e^{-i{k}_{2}t}+\mathrm{\Gamma }(t)e^{i{k}_{2}t}}{E}_{1i},$$

and subsequently the z-dependent electric field in Medium 2

(21)$${\mathbf{E}}_2(z)=\widehat{\mathbf{i}}\frac{1+\mathrm{\Gamma }(0)}{e^{-i{k}_{2}t}+\mathrm{\Gamma }(t)e^{i{k}_{2}t}}\left[e^{i{k}_2(z-t)}+\mathrm{\Gamma }(t)e^{-i{k}_2(z-t)}\right]E_0,$$

where we redefined E₀ = E_1i.

Finally, the effective electric field in Medium 2 is defined as:

(22)$${\overline{\mathbf{E}}}_{2,\text{eff}}=\frac{1}{t}{\displaystyle {\int }_0^t{\mathbf{E}}_2(z)dz}.$$

Inserting Eq. (21) in the equation above and performing the integral one obtains Eq. (7) with the notations given in Eq. (8).

Funding

European Research Council (ERC) (ERC-2014-CoG-648328); China Scholarship Council (CSC) and University College London (UCL) (201506250086).

References

1. N. Yu, P. Genevet, M. A. Kats, F. Aieta, J. P. Tetienne, F. Capasso, and Z. Gaburro, “Light propagation with phase discontinuities: generalized laws of reflection and refraction,” Science 334(6054), 333–337 (2011). [CrossRef]   [PubMed]  

2. A. V. Kildishev, A. Boltasseva, and V. M. Shalaev, “Planar photonics with metasurfaces,” Science 339(6125), 1232009 (2013). [CrossRef]   [PubMed]  

3. M. I. Shalaev, J. Sun, A. Tsukernik, A. Pandey, K. Nikolskiy, and N. M. Litchinitser, “High-efficiency all-dielectric metasurfaces for ultracompact beam manipulation in transmission mode,” Nano Lett. 15(9), 6261–6266 (2015). [CrossRef]   [PubMed]  

4. S. Jahani and Z. Jacob, “All-dielectric metamaterials,” Nat. Nanotech. 11(1), 23–36 (2016). [CrossRef]  

5. Q. Wang, E. T. Rogers, B. Gholipour, C. M. Wang, G. Yuan, J. Teng, and N. I. Zheludev, “Optically reconfigurable metasurfaces and photonic devices based on phase change materials,” Nat. Photonics 10(1), 60–65 (2016). [CrossRef]  

6. Y. Yao, R. Shankar, M. A. Kats, Y. Song, J. Kong, M. Loncar, and F. Capasso, “Electrically tunable metasurface perfect absorbers for ultrathin mid-infrared optical modulators,” Nano Lett. 14(11), 6526-6532 (2014). [CrossRef]  

7. A. Zhan, S. Colburn, R. Trivedi, T. K. Fryett, C. M. Dodson, and A. Majumdar, “Low-contrast dielectric metasurface optics,” ACS Photonics 3(2), 209-214 (2016). [CrossRef]  

8. J. Lee, M. Tymchenko, C. Argyropoulos, P.-Y. Chen, F. Lu, F. Demmerle, G. Boehm, M.-C. Amann, A. Alu, and M. A. Belkin, “Giant nonlinear response from plasmonic metasurfaces coupled to intersubband transitions,” Nature 511(7507), 65-69 (2014). [CrossRef]   [PubMed]  

9. E. Almeida, G Shalem, and Y. Prior, “Subwavelength nonlinear phase control and anomalous phase matching in plasmonic metasurfaces,” Nat. Commun. 7, 10367 (2016). [CrossRef]   [PubMed]  

10. G. Li, S. Zhang, and T. Zentgraf, “Nonlinear photonic metasurfaces,” Nat. Rev. Mater. 2, 17010 (2017) [CrossRef]  

11. A. Krasnok, M. Tymchenko, and A. Alu, “Nonlinear metasurfaces: A paradigm shift in nonlinear optics,” Mater. Today 21(1), 8–21 (2018). [CrossRef]  

12. X. Yin, Z. Ye, J. Rho, Y. Wang, and X. Zhang, “Photonic spin Hall effect at metasurfaces,” Science 339(6126), 1405–1407 (2013). [CrossRef]   [PubMed]  

13. G. Li, M. Kang, S. Chen, S. Zhang, E. Y. B. Pun, K. W. Cheah, and J. Li, “Spin-enabled plasmonic metasurfaces for manipulating orbital angular momentum of light,” Nano Lett. 13(9), 4148–4151 (2013). [CrossRef]   [PubMed]  

14. H. Wakatsuchi, S. Kim, J. J. Rushton, and D. F. Sievenpiper, “Circuit-based nonlinear metasurface absorbers for high power surface currents,” Appl. Phys. Lett. 102(21), 214103 (2013). [CrossRef]  

15. S. V. Makarov, I. S. Sinev, V. A. Milichko, F. E. Komissarenko, D. A. Zuev, E. V. Ushakova, I. S. Mukhin, Y. F. Yu, A. I. Kuznetsov, P. A. Belov, I. V. Iorsh, A. N. Poddubny, A. K. Samusev, and Yu. S. Kivshar, “Nanoscale Generation of White Light for Ultrabroadband Nanospectroscopy,” Nano Lett. 18(1), 535–539 (2018). [CrossRef]  

16. M. W. Klein, C. Enkrich, M. Wegener, and S. Linden, “Second-Harmonic Generation from Magnetic Metamaterials,” Science 313(5786), 502–504 (2006). [CrossRef]   [PubMed]  

17. W. Fan, S. Zhang, N. C. Panoiu, A. Abdenour, S. Krishna, R. M. Osgood, K. J. Malloy, and S. R. J. Brueck, “Second Harmonic Generation from a Nanopatterned Isotropic Nonlinear Material,” Nano Lett. 6, 1027 (2006). [CrossRef]  

18. S. Chen, F. Zeuner, M. Weismann, B. Reineke, G. Li, V. K. Valev, K. W. Cheah, N. C. Panoiu, T. Zentgraf, and S. Zhang, “Giant Nonlinear Optical Activity of Achiral Origin in Planar Metasurfaces with Quadratic and Cubic Nonlinearities,” Adv. Mater. 28, 2992 (2016). [CrossRef]   [PubMed]  

19. S. Liu, G. A. Keeler, J. L. Reno, M. B. Sinclair, and I. Brener, “III-V Semiconductor Nanoresonators-A New Strategy for Passive, Active, and Nonlinear All-Dielectric Metamaterials,” Adv. Opt. Mater. 10(4), 1457–1462 (2016). [CrossRef]  

20. S. Kruk, R. Camacho-Morales, L. Xu, M. Rahmani, D. A. Smirnova, L. Wang, H. H. Tan, C. Jagadish, D. N. Neshev, and Y. S. Kivshar, “Nonlinear Optical Magnetism Revealed by Second-Harmonic Generation in Nanoantennas,” Nano Lett. 17, 3914–3918 (2017). [CrossRef]   [PubMed]  

21. D. Timbrell, J. W. You, Y. S. Kivshar, and N. C. Panoiu, “A comparative analysis of surface and bulk contributions to second-harmonic generation in centrosymmetric nanoparticles,” Sci. Rep. 8, 3586 (2018). [CrossRef]   [PubMed]  

22. G. Mie, “Beiträge zur Optik trüber Medien, speziell kolloidaler Metallösungen,” Ann. Phys. (Leipzig) 25, 377 (1908). [CrossRef]  

23. A. I. Kuznetsov, A. E. Miroshnichenko, M. L. Brongersma, Y. S. Kivshar, and B. Lukyanchuk, “Optically resonant dielectric nanostructures,” Science 354, aag2472 (2016). [CrossRef]   [PubMed]  

24. P. A. Temple and C. E. Hathaway, “Multiphonon Raman spectrum of silicon,” Phys. Rev. B 7(8), 3685 (1973). [CrossRef]  

25. X. Yang and C. W. Wong, “Design of photonic band gap nanocavities for stimulated Raman amplification and lasing in monolithic silicon,” Opt. Express 13(12), 4723–4730 (2005). [CrossRef]   [PubMed]  

26. BandSOLVE™, https://www.synopsys.com.

27. A. F. Oskooi, D. Roundy, M. Ibanescu, P. Bermel, J. D. Joannopoulos, and S. G. Johnson, “MEEP: A flexible free-software package for electromagnetic simulations by the FDTD method,” Comput. Phys. Commun. 181, 687–702 (2010). [CrossRef]  

28. X. C. Zeng, D. J. Bergman, P. M. Hui, and D. Stroud, “Effective-medium theory for weakly nonlinear composites,” Phys. Rev. B 38(15), 10970–10973 (1988). [CrossRef]  

29. D. R. Smith, D. C. Vier, T. Koschny, and C. M. Soukoulis, “Electromagnetic parameter retrieval from inhomogeneous metamaterials,” Phys. Rev. E 71(3), 036617 (2005). [CrossRef]  

30. D. R. Smith and J. B. Pendry, “Homogenization of metamaterials by field averaging,” J. Opt. Soc. Am. B 23(3), 391–403 (2006). [CrossRef]  

31. A. J. Hoffman, L. Alekseyev, S. S. Howard, K. J. Franz, D. Wasserman, V. A. Podolskiy, and C. Gmachl, “Negative refraction in semiconductor metamaterials,” Nat. Mater. 6(12), 946–950 (2007). [CrossRef]   [PubMed]  

32. D. R. Smith, S. Schultz, P. Markos, and C. M. Soukoulis, “Determination of effective permittivity and permeability of metamaterials from reflection and transmission coefficients,” Phys. Rev. B 65, 195104 (2002). [CrossRef]  

33. C. G. Parazzoli, R. B. Greegor, K. Li, B. E. C. Koltenbah, and M. Tanielian, “Experimental verification and simulation of negative index of refraction using Snell’s law,” Phys. Rev. Lett. 90(10), 107401 (2003). [CrossRef]   [PubMed]  

34. D. R. Smith, W. J. Padilla, D. C. Vier, S. C. Nemat-Nasser, and S. Schultz, “Composite medium with simultaneously negative permeability and permittivity,” Phys. Rev. Lett. 84(18), 4184–4187 (2000). [CrossRef]   [PubMed]  

35. S. Roke, M. Bonn, and A. V. Petukhov, “Nonlinear optical scattering: The concept of effective susceptibility,” Phys. Rev. B 70(11), 115106 (2004). [CrossRef]  

36. A. G. F. de Beer and S. Roke, “Nonlinear Mie theory for second-harmonic and sum-frequency scattering,” Phys. Rev. B 79, 155420 (2009). [CrossRef]  

37. J. Butet and O. J. F. Martin, “Evaluation of the nonlinear response of plasmonic metasurfaces: Miller’s rule, nonlinear effective susceptibility method, and full-wave computation,” J. Opt. Soc. Am. B 33(2), A8–A15 (2016). [CrossRef]  

38. M. A. Gorlach, T. A. Voytova, M. Lapine, Y. S. Kivshar, and P. A. Belov, “Nonlocal homogenization for nonlinear metamaterials,” Phys. Rev. B 93, 165125 (2016). [CrossRef]  

39. R. W. Boyd, Nonlinear Optics (Academic, 2008).

40. B. Jalali, R. Claps, D. Dimitropoulos, and V. Raghunathan, “Light Generation, Amplification, and Wavelength Conversion via Stimulated Raman Scattering in Silicon Microstructures,” in Silicon Photonics, Top. Appl. Phys., L. Pavesi and D. J. Lockwood, (eds.) (Springer-Verlag, 2004). [CrossRef]  

41. CST Studio, https://www.cst.com.

42. J. F. McMillan, X. Yang, N. C. Panoiu, R. M. Osgood, and C. W. Wong, “Enhanced stimulated Raman scattering in slow-light photonic crystal waveguides,” Opt. Lett. 31(9), 1235–1237 (2006). [CrossRef]   [PubMed]