This is an open access article published under a Creative Commons Attribution (CC-BY) License, which permits unrestricted use, distribution and reproduction in any medium, provided the author and source are cited. Article www.acsnano.org Grating-Graphene Metamaterial as a Platform for Terahertz Nonlinear Photonics Jan-Christoph Deinert, David Alcaraz Iranzo, Raúl Pérez, Xiaoyu Jia, Hassan A. Hafez, Igor Ilyakov, Nilesh Awari, Min Chen, Mohammed Bawatna, Alexey N. Ponomaryov, Semyon Germanskiy, Mischa Bonn, Frank H.L. Koppens, Dmitry Turchinovich, Michael Gensch, Sergey Kovalev,* and Klaas-Jan Tielrooij* Cite This: ACS Nano 2021, 15, 1145−1154 Read Online ACCESS Metrics & More Article Recommendations * sı Supporting Information ABSTRACT: Nonlinear optics is an increasingly important field for scientific and technological applications, owing to its relevance and potential for optical and optoelectronic technologies. Currently, there is an active search for suitable nonlinear material systems with efficient conversion and a small material footprint. Ideally, the material system should allow for chip integration and room-temperature operation. Two-dimen- sional materials are highly interesting in this regard. Particularly promising is graphene, which has demonstrated an exceptionally large nonlinearity in the terahertz regime. Yet, the light−matter interaction length in two-dimensional materials is inherently minimal, thus limiting the overall nonlinear optical conversion efficiency. Here, we overcome this challenge using a metamaterial platform that combines graphene with a photonic grating structure providing field enhancement. We measure terahertz third-harmonic generation in this metamaterial and obtain an effective third-order nonlinear susceptibility with a magnitude as large as 3 × 10−8 m2/V2, or 21 esu, for a fundamental frequency of 0.7 THz. This nonlinearity is 50 times larger than what we obtain for graphene without grating. Such an enhancement corresponds to a third-harmonic signal with an intensity that is 3 orders of magnitude larger due to the grating. Moreover, we demonstrate a field conversion efficiency for the third harmonic of up to ∼1% using a moderate field strength of ∼30 kV/cm. Finally, we show that harmonics beyond the third are enhanced even more strongly, allowing us to observe signatures of up to the ninth harmonic. Grating-graphene metamaterials thus constitute an outstanding platform for commercially viable, CMOS-compatible, room-temperature, chip- integrated, THz nonlinear conversion applications. KEYWORDS: terahertz, harmonics, graphene, nonlinear, field enhancement, metamaterial T he ability to perform nonlinear optical conversion has relevance for a broad range of current and future technologies, including ultrashort pulse generation, advanced spectroscopy, optical information processing and made in recent years, with suitable material systems operating in the visible and near-infared parts of the electromagnetic spectrum, for example, using indium tin oxide with a thickness of tens of nanometers.5,6 However, currently, no commercially viable nonlinear material system is available for incident light storage, telecommunications, light harvesting, bioimaging, integrated optics, and quantum technologies.1−4 For many of in the terahertz region (roughly 0.3−30 THz) of the these applications, it is crucial to find nonlinear optical material electromagnetic spectrum. This technological gap needs to systems that allow for integration into existing nanophotonic be filled because the THz region holds great promise for both platforms and for operation with low power consumption. Thus, ideal nonlinear optical materials should be compatible Received: September 25, 2020 with complementary metal−oxide−semiconductor (CMOS) Accepted: November 25, 2020 technology, have a small material footprint, and have a large Published: December 11, 2020 nonlinear susceptibility, such that large nonlinear effects occur for small incident optical fields. Furthermore, they should ideally operate at room temperature. Steady progress has been © 2020 American Chemical Society https://dx.doi.org/10.1021/acsnano.0c08106 1145 ACS Nano 2021, 15, 1145−1154 ACS Nano www.acsnano.org Article Figure 1. Enhanced third-harmonic generation in grating-graphene metamaterial. (a) Schematic representation of the measurement configuration, with a multicycle THz waveform with fundamental frequency f incident on a grating-graphene metamaterial sample (top) or a bare graphene sample (bottom). In the metamaterial, field enhancement occurs inside the metal gap, which leads to nonlinearly enhanced generation of third-harmonic signal. (b,c) Measured THz field strength in the time (b) and frequency (c) domain for the grating-graphene metamaterial, with an incident field strength of 12 kV/cm. (d,e) Measured THz field strength in the time (d) and frequency (e) domain for the bare graphene sample with an incident field strength of 13.6 kV/cm. Clearly, a significantly more intense harmonic signal is created in the graphene metamaterial sample. scientific and technological applications, such as imaging, graphene was found to be especially large: |χ(3)| ≈ 10−9 m2/ industrial quality control, gas sensing, and wireless communi- V2.16 This value is many orders of magnitude larger than its cations. All of the above would benefit tremendously from the χ(3) in the near-infrared range, which is below ∼10−16 m2/V2.21 availability of integrated nonlinear optical devices in the THz The main reason for the strong THz nonlinearity is that it is range. the result of a different nonlinearity mechanism than the one in It is particularly appealing to use nonlinear optical materials the near-infrared. Namely, the nonlinearity originates from the that are thinner than the wavelength of light. In this case, collective thermodynamic response of the quasi-free carriers in different optical waves propagating through the material graphene. Their response to an incident THz pulse can be automatically remain in phase, thus removing the need for outlined as follows (see also refs 16, 22, and 23): as a first step, complex engineering to achieve phase matching. As a result, graphene conduction-band electrons absorb the energy from nonlinear optical conversion in ultimately thin, two-dimen- incident THz light through Drude absorption. This energy is sional materials has received ample attention, as reviewed then distributed among all electrons, leading to a thermal recently in ref 7. However, the challenge associated with using quasi-equilibrium within the electronic system, as it has not yet subwavelength materials is thateven with extraordinary large thermalized with the lattice. This is because the THz nonlinear coefficientsthere is simply not much material to oscillation period and the electron cooling time are on the interact with, thus usually leading to low nonlinear conversion order of a picosecond and thus much longer than the electron efficiencies. A promising strategy to overcome this challenge thermalization time, which is below 100 fs. The electronic had already been identified before the advent of two- temperature increase of the quasi-equilibrium Fermi−Dirac dimensional materials. Namely, photonic structures can be distribution is substantial due to the small electronic heat used to enhance light−matter interactions and thereby increase capacity of graphene.24 Since the THz conductivityand the nonlinear conversion efficiency. Typical photonic struc- tures that have been used include nanoparticles,8 hole therefore THz absorptiondecreases for increasing electron arrays,9,10 nanocavity gratings,11 nanocavity−antenna sys- temperature (and Fermi energy above ∼0.1 eV),22,25,26 this tems,12 split-ring resonators,13,14 and nanoantennas.4,6,15 gives rise to a large, ultrafast thermodynamic THz nonlinearity. Combining a photonic structure with a nonlinear material Quasi-monochromatic THz excitation leads to ultrafast has, for example, made it possible to obtain a surface-enhanced modulation of the graphene THz conductivity, and hence of Raman spectrum of a single molecule8 and to observe third- the driven THz currents, within each cycle of the driving THz harmonic generation from a single semiconductor nano- field due to the electronic heating−cooling dynamics. This, in particle.4 turn, leads to a nonlinear response that manifests as harmonic Nonlinear effects in two-dimensional materials are partic- generation at odd-order multiples of the frequency of the ularly strong for light in the THz region. Indeed, two- driving field. Given this strong nonlinearity and the fact that dimensional, gapless Dirac materials, where electrons obey a graphene can be integrated with CMOS-based chips, as linear energy-momentum dispersion relation, have shown very demonstrated in ref 27, it is a highly promising candidate large nonlinear coefficients in the THz range. In particular, material for commercially viable, chip-integrated, THz non- THz third-harmonic generation has been observed using linear conversion applications. One important hurdle, however, monolayer graphene,16 using the Dirac surface states in the that needs to be overcome is to ensure sufficient light−matter topological insulator Bi2Se3,17,18 and using the Dirac semimetal interaction between incident THz light and a monolayer thin Cd3As2.19,20 The THz third-order nonlinear susceptibility of material. 1146 https://dx.doi.org/10.1021/acsnano.0c08106 ACS Nano 2021, 15, 1145−1154 ACS Nano www.acsnano.org Article Figure 2. Nonlinearity of grating-graphene metamaterial vs bare graphene. (a) Comparison of third-harmonic intensity for grating-graphene metamaterial sample A (red squares) and bare graphene sample B (blue circles), as a function of peak field strength of the incident THz light. The blue dashed line is the result of split-step simulations based on thermodynamic nonlinearity without field enhancement. The red dashed line represents the same simulation, now with field-enhancement factor M = 5 and active area (1 − η) = 10%. (b) Extracted third- (3) order nonlinear susceptibility χeff as a function of incident peak field strength for grating-graphene metamaterial sample A (red squares) and bare graphene sample B (blue circles). The dashed lines show the results of the thermodynamic model. The decreasing trends with incident field strength indicate the occurrence of saturation effects. For low field strengths, the metal grating leads to a factor 40−50 increase in χeff (3) , corresponding to more than 3 orders of magnitude enhanced harmonic intensity. The red and blue shaded regions in both panels indicate the onset of saturation effects not captured by the model. In this work, we demonstrate a metamaterial system order to be able to resolve this signal well within the dynamic consisting of graphene and a metallic grating, where the range of the detectors. Clearly, the measured oscillations are grating provides field enhancement, thus leading to a strong faster than those of the incident light. The Fourier transform increase of the THz nonlinear light−matter interaction. We (see Methods and Supporting Information) in Figure 1c shows examine the generation of odd THz harmonics and observe that these fast oscillations correspond to the generated third that the third-harmonic intensity can be enhanced by more harmonic at 3f = 2.1 THz. than 3 orders of magnitude due to the integration of the metal We first compare the nonlinear conversion efficiency of two grating. We discuss quantitatively how this enhancement is samples, both consisting of wet-transferred monolayer achieved and what the limitations of this approach are and the graphene on quartz with an estimated Fermi energy of ∼0.23 perspectives for further enhancing THz nonlinearities with this eV (see Methods and Supporting Information for details). photonic approach. Finally, we show a field conversion Sample A has a gold metallic grating on top of the graphene, efficiency of 1% and signatures of higher harmonics, including separated by a 2 nm thick layer of Al2O3. The metal grating the ninth harmonic, using a moderate incident field strength of contains metallic stripes with a metal width of wmetal = 18 μm, ∼30 kV/cm. separated by gaps with a width of wgap = 2 μm. This w corresponds to a duty cycle η = w metal of 90%. Thus, 90% RESULTS AND DISCUSSION +w metal gap Grating-Enhanced THz Third-Order Nonlinearity. A of the graphene area is covered by metal, and there is 10% natural and proven method to study the nonlinear properties active graphene area. This duty cycle is the parameter that of a material is to quantify the intensity or field strength of governs the field enhancement and absorption in the grating- generated harmonic signal compared to the intensity or field graphene metamaterial samples, as we will show. Sample B also strength of incident fundamental light. This method has contains 2 nm of Al2O3 on top of graphene, while not already been applied successfully to graphene,16 and therefore, containing a metal grating. As a first indication that our we follow the same experimental procedure. Briefly (see approach of using a grating-graphene metamaterial to enhance Methods and Supporting Information for details), this method nonlinear conversion works, we use a small incident field consists of focusing a narrow-band THz waveform, generated strength of Ef = 12 kV/cm and observe clear third-harmonic by the accelerator-based TELBE beam facility,28 onto a sample signal for the grating-graphene metamaterial sample A (Figure and measuring the electric field of the transmitted waveform as 1b,c). However, for the bare graphene sample B with a slightly a function of time through electro-optic sampling. In Figure 1a, larger incident field strength of 14 kV/cm, there is no clearly we show a measured incident THz waveform with fundamental observable third-harmonic signal in the time trace (see Figure frequency f = 0.7 THz in the time domain, together with a 1d). In the Fourier spectrum (see Figure 1e), we observe a schematic of a grating-graphene metamaterial sample and a small signal just above the noise floor. Clearly, the bare graphene sample. Figure 1b shows the signal in the time metamaterial sample has a much stronger nonlinearity, even domain after passing through the metamaterial sample and though only 10% of the graphene is part of the active area. We spectral filters (two band-pass filters at 2.1 THz). These band- note that the metallic grating itself does not absorb much THz pass filters suppress spectral components at the fundamental light and does not generate any significant third-harmonic frequency, while transmitting the third-harmonic signal, in signal (see Supporting Information). 1147 https://dx.doi.org/10.1021/acsnano.0c08106 ACS Nano 2021, 15, 1145−1154 ACS Nano www.acsnano.org Article In order to quantify the nonlinear properties of our grating- for the simulation. This produces the blue dashed lines in graphene metamaterial, we compare how the field strength of Figure 2, which describe the experimental data quite well up to the third-harmonic signal E3f scales with the incident field an incident field strength of ∼20 kV/cm. For higher field strength Ef (see Figure 2a). These data points are obtained by strengths, there is more saturation in the experimental data taking the Fourier transform of transmitted time traces as in than what is predicted by the model. We will discuss saturation Figure 1 and taking the peak of the power spectrum around the effects in more detail in the next section. For the simulations of frequency of the third-harmonic peak at 2.1 THz. We observe metamaterial sample A, we use the exact same input that, for the lowest field strengths, E3f is more than an order of parameters as for sample B, with the difference that we use magnitude higher in the metamaterial sample than in bare as incident field strength M × Ef, where M is the field- graphene, whereas the third-harmonic signals are more similar enhancement factor. We estimate an average field-enhance- for the two samples at the highest field strengths. We extract ment factor M of ∼5 inside the gap of our grating, using the third-order susceptibility using16 rigorous coupled wave analysis (RCWA), following ref 34 (see Supporting Information). The resulting third-harmonic signal (3) E3f 4c·n3f would correspond to a sample where the whole graphene area |χeff |= Ef3 3πfd (1) experiences 5-fold field enhancement, whereas in our case, the active area is only 1 − η = 10%. Therefore, we multiply the where c is the speed of light and d = 0.3 nm is the thickness of simulation results with a factor 1 − η , where the square root graphene. For the refractive index of graphene at the third- harmonic frequency, we conservatively use n3f = 10, as in ref comes from the fact that we consider field rather than intensity. 16. This value can be up to an order of magnitude larger, which The results are the red dashed lines in Figure 2, which match would make the nonlinear coefficient an order of magnitude the data for the lowest field strengths. Thus, the experimental larger. The results in Figure 2b show that for low incident field results are in agreement with a grating-induced field enhance- strength the grating-graphene metamaterial reaches an ment of a factor ∼5. Importantly, the third-harmonic signal is increased in a (3) experimental value of |χeff | = 3 × 10−8 m2/V2 (using n3f = nonlinear fashion by the field enhancement. Indeed, we find 10). In esu units, this corresponds to29 that for the lowest incident field strengths, the model predicts (3) |χeff |[esu] = 9 × 109 (3) |χeff | = 21 esu. Expressed as a surface an increase of the third-harmonic field strength by a factor 4π ∼40, corresponding to an intensity enhancement factor of nonlinearity, we obtain |χ2D (3) (3) | = |χeff |·d = 10−17 m3/V2. We >1000. We can understand these numbers analytically, as for can also express the nonlinearity as the nonlinear refractive the nth harmonic, we expect index n2, where the intensity-dependent refractive index is given as1 n = n0 + n2 × ITHz. Here n0 is the linear refractive Enf(grating) index and ITHz is the intensity of the incoming THz light. We = M n· 1 − η Enf(nograting) (2) estimate an n2 well above 100 cm2/GW (see Methods). Finally, we note that the sign of the THz nonlinearity, and (3) For n = 3, M = 5, and η = 90%, the analytical model gives the therefore of χeff and n2, is negative, as the refractive index same factor of ∼40 as the numerical model based on the decreases upon THz excitation (see also ref 16), in agreement thermodynamic nonlinearity. This factor is in good agreement with the sign of the nonlinearity in the near-infrared region.30 (3) with the experimental observation of a |χeff | that is ∼50 times Comparing the nonlinear susceptibility of the grating- larger for the grating-graphene metamaterial sample A, in graphene metamaterial sample A with graphene without comparison with the bare graphene sample B. The small grating (sample B), we find a value that is ∼50 times larger. (3) discrepancy could come from a larger field-enhancement factor The obtained value of |χeff | is also more than 3 orders of M or from the fact that neither the thermodynamic nor the magnitude larger than the highest THz nonlinearities reported analytical model take into account that the grating could also for quantum well systems.31−33 We can also compare the THz lead to more efficient out-coupling of harmonic light to the far- nonlinearity of our grating-graphene metamaterial with the field, where it is detected. We note that in the graphene-filled near-infared nonlinearities found very recently for thin films of gap region, and in the perturbative low-field regime, the third- indium tin oxide in a special frequency region with epsilon near harmonic intensity is enhanced by a factor M6 (see eq 1), zero, with a nonlinear refractive index n2 of ∼0.1 cm2/GW which is more than 4 orders of magnitude. without field enhancement5 and ∼4 cm2/GW with nano- Saturation of Harmonic Generation. We now discuss antenna-induced enhancement.6 These nonlinear refractive the occurrence of saturation effects, which follow clearly from indices are lower than what we find for our grating-graphene Figure 2b, where without saturation one would obtain a flat metamaterial in the THz regime. Overall, we can conclude that line. We first discuss two saturation mechanisms that are the THz nonlinearity of our grating-graphene metamaterial is captured by the thermodynamic model. The first reason for large, compared to any nonlinear optical (meta)material at any saturation is the occurrence of saturable absorption: if the wavelength, operating at room temperature. carrier temperature increases, the graphene conductivity We can understand the enhancement of the nonlinearity decreases and the absorption of THz light is reduced (see induced by the metal grating quantitatively by incorporating also Supporting Information). Eventually, the carrier temper- field enhancement into the thermodynamic model of refs 16 ature can become so high that virtually no incident THz light is and 22. For bare graphene sample B, we straightforwardly use absorbed, as shown in ref 22. The second reason that the incident field strength Ef, the modal carrier density saturation occurs is thatwith constant absorptionthe obtained from Raman measurements (see Supporting In- increase in carrier temperature of graphene scales sublinearly formation), an estimated momentum scattering time, and with the intensity of the incident light, which is the result of standard picosecond cooling dynamics as the input parameters the electronic heat capacity scaling with temperature.35 This 1148 https://dx.doi.org/10.1021/acsnano.0c08106 ACS Nano 2021, 15, 1145−1154 ACS Nano www.acsnano.org Article Figure 3. Effect of grating orientation and duty cycle. (a) Field strength of the created third-harmonic signal as a function of orientation of the metal grating with respect to the polarization of the incoming fundamental THz light (red squares; sample A; incident field strength of ∼32 kV/cm), as illustrated in the schematic insets. The dashed line is a sin2 function to guide the eye. (b) Schematic representation (top) and optical images (bottom) of the four different grating-graphene areas of sample C, with duty cycles of 24, 44, 66, and 79%. Metal and gap widths are to scale. (c) Field strength of the created third-harmonic signal as a function of duty cycle using sample C (red squares; incident peak field strength of ∼14 kV/cm), showing an increase, followed by saturation, for increasing duty cycle. The thermodynamic model (dashed line) matches the experimental results for the lowest duty cycles before saturation effects occur, as indicated by the red shaded region. The inset shows the absorption of fundamental light vs duty cycle, showing larger absorption at higher duty cycle. This is obtained using (1 − transmitted fundamental intensity) and multiplying this by a constant, such that the area with a duty cycle of 24% has a (1 − transmission) ≈ absorption of 33%, taken from RCWA simulations. (d) RCWA simulations of the absorption of grating-graphene structures as a function of Fermi energy and duty cycle, predicting higher absorption and thus higher field-enhancement factor M for higher duty cycle. The increasing trend with Fermi energy at moderate duty cycles comes from the intrinsically larger THz absorption when there are more graphene charges to be accelerated. sublinear behavior starts to play a role when the carrier Supporting Information). With heat accumulating in both the temperature increase ΔT is larger than the initial carrier, or electronic and phononic systems, phonon bottleneck effects lattice, temperature (300 K). We note that it is the increase in lead to a longer cooling time for the hot carrier system. We carrier temperature that gives rise to the THz nonlinearity.16,22 have simulated the effect of slower cooling and find that a Thus, the sublinear scaling of carrier temperature with light factor 10 increase in cooling time leads to almost a factor 10 intensity gives rise to a less strong nonlinearity. decrease in harmonic generation (see Supporting Informa- There are also saturation effects that are not captured by our tion). Furthermore, it gives rise to a red shift of the harmonic thermodynamic model, thus leading to a discrepancy between signal, which we indeed observe (see Supporting Information). experimental data points and the results of the thermodynamic At these elevated carrier temperatures, increased phonon model (see the blue and red shaded areas in Figure 2). We first scattering can furthermore decrease the carrier mobility, thus discard the occurrence of light-induced damage to the grating giving rise to even stronger saturable absorption. Therefore, we or the graphene, as we did not observe any nonreversible ascribe the strong saturation effects we observe in the grating- effects of the generated harmonic signal when increasing and graphene metamaterial, even for field strengths on the order of decreasing incident field strength, and confirmed that there 10 kV/cm, to the comparably high carrier temperatures that we was no damage by inspection with an optical microscope. We reach due to field enhancement. note that the damage observed in ref 13 inside their split-ring Photonic Enhancement. We proceed with exploring how resonator gap occurred for a peak field strength of ∼4 MV/cm, the intensity of the generated third-harmonic signal depends an order of magnitude larger than our highest local field inside on the photonic enhancement of the metal grating and what the grating gap, which was M × 60 = 300 kV/cm. Having are the design considerations to optimize the grating. First, we discarded any damage, we consider saturation related to the verify the proper operation of the grating by varying its carrier temperature of graphene that goes beyond the orientation with respect to the polarization of the incident thermodynamic model. We point out that for the metamaterial field. The measurements in Figure 1 were performed with the sample the carrier temperature reaches a value above 2000 K metal stripes perpendicular to the incoming polarization, where already for an incident field strength of ∼20 kV/cm (see enhancement of field inside the gap of the grating structure is 1149 https://dx.doi.org/10.1021/acsnano.0c08106 ACS Nano 2021, 15, 1145−1154 ACS Nano www.acsnano.org Article expected. Figure 3a shows that by rotating the sample such that Fermi energy of graphene, where we note that a smaller Fermi the orientation of the grating with respect to the incoming energy means less intrinsic THz absorption. Our simulations THz polarization changes over 90°, the amount of generated indicate that the combination of graphene with a Fermi energy harmonic signal is reduced to a value below our experimental approaching the Dirac point and a metal grating with a duty noise floor. This is also as expected because, for parallel cycle approaching 100% still leads to absorption close to 50%. orientation of the polarization with respect to the metal stripes, As an outlook of what is possible with our nonlinear grating- no field enhancement occurs. Moreover, the grating in this case graphene metamaterial platform, we calculate the expected behaves as a THz filter that reduces THz transmission to the nonlinearity for a grating with a duty cycle of 99% and graphene film behind it. graphene with a Fermi energy of 0.1 eV. Such a duty cycle can In order to study how field enhancement depends on the be achieved using a gap size of wgap = 200 nm and (wmetal + duty cycle of the grating, we use a grating-graphene wgap) = 20 μm. This is possible using nanofabrication metamaterial sample (sample C), which contains four areas techniques, such as electron beam lithography or focused ion with metal gratings with different duty cycles, varying from 24 beam milling. We estimate that we can reach M = 30 for η = to 79% (see Figure 3b). A larger duty cycle corresponds to a 99% and EF = 0.1 eV. This would increase the third-harmonic smaller gap in between the metal stripes, which should lead to intensity by a factor M6 × (1 − η) ≈ 107 and increase the larger field enhancement. Figure 3c shows that, with increasing third-order THz nonlinearity to above 10−6 m2/V2. duty cycle, the amount of generated harmonic indeed Overall Conversion Efficiency and Ninth-Harmonic increases, although it saturates for the largest duty cycles. Generation. In addition to optimizing the nonlinear The reason that we do not see any further increase in detected susceptibility, as we just described, an important goal toward third harmonic above η = 70% (see Figure 3c) is the applications is to obtain a large nonlinear conversion efficiency occurrence of saturation effects, as discussed above. In order to for relatively small incident field strength. We therefore test gain further insights, we have performed RCWA simulations as sample A (with 90% duty cycle) using an experimental a function of duty cycle and Fermi energy of graphene (see geometry where we directly obtain the conversion efficiency. Figure 3d). The results show that the total amount of THz To this end, we use a table-top THz source, producing a absorption increases with increasing duty cycle before moderate peak field strength of ∼30 kV/cm. Briefly (see decreasing quickly upon approaching η = 100%. This total Methods for details), we use an amplified Ti:sapphire laser absorption is the sum of the absorption of the active graphene system incident on a LiNbO3 nonlinear crystal, together with area inside the gap of the grating and of the inactive area that is two 0.3 THz band-pass filters, to create a multicycle THz covered by metal. The field in the graphene region below the waveform in a narrow band around 0.3 THz. As a detection metal stripes is very small, so absorption there is negligible (see crystal, we use a 100 μm thick ZnTe crystal with ⟨110⟩ cut on Supporting Information). For our samples with Fermi energy top of 2 mm ⟨100⟩ ZnTe.36 This gives us a flat spectral of ∼0.23 eV, we expect maximum absorption for a duty cycle response and eliminates reflections, meaning that the of 80−90%. This behavior can be understood from the trade- sensitivities for the fundamental at 0.3 THz and the third off between increasing field enhancement vs decreasing active harmonic at 0.9 THz are identical. Crucially, we measure the graphene area with increasing duty cycle. We measured the transmitted THz waveform without any filters behind the transmitted fundamental intensity for sample C with varying sample, thus transmitting both fundamental and created duty cycle and used it to examine how the absorption ≈ (1 − harmonic signal. The Fourier transform is shown in Figure transmitted fundamental) scales with duty cycle (see inset of 4a and clearly contains both fundamental fields at 0.3 THz and Figure 3c). We confirm that the overall THz absorption (active third harmonic at 0.9 THz. We directly obtain the ratio plus inactive area) increases moderately with increasing duty between generated third-harmonic and transmitted fundamen- cycle. As the active area decreases with increasing duty cycle, tal fields and find a value of 1.6%, whichafter correcting for this indicates increasingly strong field enhancement inside the an estimated ∼45% absorption of THz fundamental in the gap. metamaterialcorresponds to a nonlinear conversion effi- The absorption simulations suggest that a grating-graphene ciency of ∼1%. This is a highly encouraging value for the field metamaterial, despite the small active area, can reach almost conversion efficiency of THz third-harmonic generation using 50% absorption of all incident THz photons onto the such a modest field strength (∼30 kV/cm). metamaterial for monolayer graphene with a standard mobility As a final experiment, we examine harmonics beyond the of 2000 cm2/Vs. This is ∼3 times more than the amount of third harmonic. As grating-induced field enhancement absorption in bare graphene with a Fermi energy of 0.23 eV, as increases harmonic generation in a nonlinear way (see eq 1), estimated for sample A (see Supporting Information). With a this effect will be even more pronounced for higher harmonics. duty cycle of η = 90%, that is, 10% active area, the simulations Analytically, we would expect the intensity of the n = ninth give an overall absorption of 45% in the metamaterial sample, harmonic, for example, to be enhanced by 11 orders of which corresponds to an effective absorption of 450% in the magnitude using our grating with a duty cycle of 90% and field- active graphene area. This is a photonic absorption enhance- enhancement factor M = 5. This, however, does not take any ment of ∼26, which corresponds to a field enhancement of M saturation effects into consideration. Therefore, we obtain a ∼ 5 inside the gap, which is in good agreement with the more accurate estimate using the thermodynamic split-step electric field simulations in the Supporting Information. For calculations of THz harmonic generation following ref 16. We the lowest duty cycle, we studied experimentally, η = 24%, the use a modest incident field strength of 30 kV/cm in order to overall THz absorption is still around 30%. However, the field avoid the occurrence of saturation effects that are not captured enhancement M inside the gap is only ∼1.5. Thus, in order to by our model. The simulation results in Figure 4b predict that maximize the nonlinear susceptibility, the grating-graphene the ninth-harmonic intensity for our grating-graphene non- metamaterial would have a duty cycle approachingwithout linear metamaterial will be enhanced by more than 7 orders of reaching100%. The second parameter to optimize is the magnitude compared to that of bare graphene. To verify this 1150 https://dx.doi.org/10.1021/acsnano.0c08106 ACS Nano 2021, 15, 1145−1154 ACS Nano www.acsnano.org Article already large THz nonlinearity of graphene. For a grating with a duty cycle of 90%, we have obtained an effective nonlinear (3) susceptibility χeff as high as −3 × 10−8 m2/V2 and a third- harmonic generation conversion efficiency of ∼1% using a moderate incident field strength of ∼30 kV/cm. We qualitatively understand the occurrence of this large non- linearity as being the result of the incident field being enhanced inside the gaps of the grating. Due to the nonlinear nature of harmonic generation, the intensity of the third-harmonic signal is increased by a factor much larger than the field enhance- ment. Quantitatively, we find agreement using a field- enhancement factor of 5, as obtained from RCWA calculations, and calculating harmonic generation analytically in the perturbative regime and numerically using a thermodyamic model of the THz nonlinearity beyond the perturbative regime. We have shown that, for incident field strengths above ∼10 kV/cm, saturation effects related to the high carrier temperature occur. Low incident field strengths are in any case most relevant, as commercial applications require low power consumption. Based on our understanding of the metal grating, we have shown perspectives for further improving the nonlinear conversion efficiency of the grating-graphene system, in particular, by increasing the duty cycle and decreasing the Fermi energy. Finally, we have shown that the nonlinear enhancement allows for creating higher harmonics (up to the Figure 4. Conversion efficiency and ninth harmonic. (a) Spectral ninth). The grating-graphene metamaterial that we have content of the THz waveform after grating-graphene metamaterial introduced here relies on standard, CMOS-compatible sample A using a table-top source with 30 kV/cm at 0.3 THz. We fabrication technology using chemical vapor deposition did not use any spectral filters in the detection pulse, thus allowing (CVD)-grown graphene and has a very small material footprint the amount of transmitted fundamental fields and generated harmonic in the spectrum to be quantitatively observed, giving a and sufficiently high nonlinear conversion efficiency to ensure field conversion efficiency of ∼1%. (b) Split-step simulation results low power consumption. We thus conclude that it is an of the effect of a grating-induced enhancement factor M = 5 on the outstanding platform for commercially viable, chip-integrated, possibility to observe higher harmonics. The incident field is 30 THz nonlinear conversion applications. kV/cm at 0.3 THz. The blue (red) line is without (with) field enhancement. At the fundamental frequency, the intensity is METHODS enhanced by M2 = 25, whereas, due to the higher-order nonlinear processes that give rise to higher harmonics, field confinement Sample Fabrication. All samples were prepared starting with 1 leads to even larger enhancement of harmonic intensity, for cm × 1 cm Infrasil quartz substrates with excellent THz transparency. example, around 7 orders of magnitude for the ninth harmonic. (c) For samples A and B, we then used standard poly(methyl Experimental results using our table-top source with incident field methacrylate) (PMMA)-assisted wet transfer of CVD-grown strength of ∼30 kV/cm and grating-graphene metamaterial sample graphene (Grolltex) of roughly 1 cm × 1 cm. This was followed by A. We put high-pass filters in the detection path in order to focus atomic layer deposition (ALD) of 2 nm of Al2O3. For sample A, we on higher harmonics, thus suppressing strongly the fundamental then used optical lithography and thermal evaporation of 50 nm of and third harmonic. We observe very clear signatures of the fifth gold with titanium as the sticking layer. Samples A and B, which are and seventh harmonic and even the ninth harmonic. The gray used for comparing the grating-induced harmonic enhancement, were shaded area indicates the estimated experimental noise floor. prepared concurrently, under identical conditions with the exact same source materials, thus ensuring that their key propertiesFermi energy and mobilityare nearly identical. For sample C, we first experimentally, we use our table-top setup, now with two high- fabricated the metal grating, followed by ALD of 2 nm of Al2O3 and pass filters behind sample A and a GaP detection crystal, which then standard PMMA-assisted wet transfer of CVD-grown graphene. allows for a larger acceptance THz spectral bandwidth, as The grating areas had a typical size of about 5 mm × 5 mm, that is, significantly larger than the spot size. See Supporting Information for compared to the previously used ZnTe crystal. Figure 4c shows a schematic overview of the samples and Raman characterization of the obtained spectrum with very clear fifth and seventh sample B. We note that RCWA calculations show no significant harmonic. There is even a discernible signature of the ninth difference between samples where the graphene is placed on top or harmonic just above our estimated noise floor. We remark that below the metal grating (see Supporting Information). the observation of ninth harmonic signal with incident light in Experimental Setups. Measurements were all performed at the the near-infrared in the strongly nonlinear epsilon near zero TELBE facility at Helmholtz-Zentrum Dresden-Rossendorf, at room material CdO required a much larger field strength, namely, >1 temperature and ambient conditions.28 See Supporting Information MV/cm.37 This highlights the exceptional strength of the THz for a schematic of the experimental setup. The super-radiant source nonlinearity of our material system. was tuned to provide multicycle, monochromatic, linearly polarized THz pulses with a frequency of 0.7 THz (20% fwhm bandwidth) at a repetition rate of 100 kHz. The field strength was tuned using a pair CONCLUSION of wiregrid polarizers. The incident THz light was focused down to a We have demonstrated a grating-graphene metamaterial, where spot size with a fwhm of 870 μm. We used two band-pass filters at 2.1 grating-induced field enhancement significantly increases the THz behind the sample to reduce the transmitted fundamental field. 1151 https://dx.doi.org/10.1021/acsnano.0c08106 ACS Nano 2021, 15, 1145−1154 ACS Nano www.acsnano.org Article The transmitted THz waveform, containing mainly re-emitted third- ASSOCIATED CONTENT harmonic signal, was detected with electro-optic sampling (EOS) in a 1.9 mm thick ZnTe nonlinear crystal using a synchronized Ti:sapphire * sı Supporting Information laser system with 100 fs pulses.38 In the analysis of the transmitted The Supporting Information is available free of charge at signals, we took into account the properties of the two filters, as well https://pubs.acs.org/doi/10.1021/acsnano.0c08106. as the wavelength-dependent EOS response of ZnTe. For more Figures of schematic of experimental THz harmonics details, see ref 16. setups, data analysis with window function, sample For the measurements in Figure 4, we used a Ti:sapphire amplified description and characterization, grating-graphene meta- laser system to pump a nonlinear LiNbO3 crystal under a tilted pulse material vs bare metal grating, RCWA simulations of front configuration39 to generate THz radiation with a broadband, field enhancement, RCWA simulations of absorption linearly polarized THz waveform at a repetition rate of 1 kHz. We enhancement, saturation effects, effect of cooling then used a band-pass filter at 0.3 THz to create a spectrally more dynamics, fundamental and harmonic signals in narrow, multicycle THz waveform. We detected the transmitted THz radiation in Figure 4a, without passing through any filters, using a 100 frequency domain (PDF) μm thick ZnTe crystal with ⟨110⟩ cut on top of 2 mm ⟨100⟩ ZnTe.36 For the data in Figure 4c, we performed EOS using a 300 μm thick AUTHOR INFORMATION GaP crystal, and we added two high-pass filters with cutoff frequencies Corresponding Authors at 1.4 and 2.2 THz in order to resolve higher harmonics. Split-Step Time-Domain Simulation of Thermodynamic THz Sergey Kovalev − Helmholtz-Zentrum Dresden-Rossendorf, Nonlinearity. We treat the nonlinear time-dependent response of Dresden 01328, Germany; Email: s.kovalev@hzdr.de graphene interacting with an incident THz field as described in detail Klaas-Jan Tielrooij − Catalan Institute of Nanoscience and in refs 16, 22, and 23. Briefly, the split-step approach consists of Nanotechnology (ICN2), BIST & CSIC, Bellaterra calculating the response of the system in very small temporal steps, in (Barcelona) 08193, Spain; orcid.org/0000-0002-0055- our case, steps of 1 fs. After the arrival and absorption of the first part 6231; Email: klaas.tielrooij@icn2.cat of the THz pulse by graphene, the amount of heat contained in the electronic system increases. This corresponds to an increased electron Authors temperature, which is associated with a lower THz conductivity, that Jan-Christoph Deinert − Helmholtz-Zentrum Dresden- is, lower THz absorption. This therefore affects the interaction Rossendorf, Dresden 01328, Germany; orcid.org/0000- between graphene and the next part of the incoming THz pulse. At 0001-6211-0158 the same time, the electronic system also starts to cool down, thus David Alcaraz Iranzo − ICFO - Institut de Ciències increasing again the THz absorption. These THz-induced heating− Fotòniques, The Barcelona Institute of Science and cooling cycles of the electronic system, which lead to a temporal Technology, Castelldefels (Barcelona) 08860, Spain modulation of the THz absorption, are the basis of the Raúl Pérez − Catalan Institute of Nanoscience and thermodynamic THz nonlinearity. In our split-step simulation, we Nanotechnology (ICN2), BIST & CSIC, Bellaterra repeat the calculation of the system response every femtosecond in a (Barcelona) 08193, Spain looped numerical fashion. This allows us to obtain the buildup of the Xiaoyu Jia − Max-Planck-Institut für Polymerforschung, heat function in time and extract the transmitted field until the whole Mainz 55128, Germany driving THz field is propagated through the graphene layer. Finally, Hassan A. Hafez − Fakultät für Physik, Universität Bielefeld, we use Fourier transformation to obtain the fundamental and Bielefeld 33615, Germany harmonic components of the transmitted field. We use cooling dynamics with characteristic time scales in the picosecond regime, as Igor Ilyakov − Helmholtz-Zentrum Dresden-Rossendorf, in refs 16 and 23. We furthermore apply a phenomenological Dresden 01328, Germany parametrization of the effect of the amount of absorbed heat in the Nilesh Awari − Helmholtz-Zentrum Dresden-Rossendorf, electronic system on the THz conductivity (and thus absorption), as Dresden 01328, Germany in refs 16 and 23. Min Chen − Helmholtz-Zentrum Dresden-Rossendorf, Dresden Calculation of Nonlinear Refractive Index. For Kerr-like 01328, Germany optical nonlinearities, the nonlinear refractive index is given by1 Mohammed Bawatna − Helmholtz-Zentrum Dresden- Rossendorf, Dresden 01328, Germany (3) 3χeff Alexey N. Ponomaryov − Helmholtz-Zentrum Dresden- n2 = for |n2|ITHz /n0 ≪ 1 Rossendorf, Dresden 01328, Germany 4n3f2 ϵ0c (3) Semyon Germanskiy − Helmholtz-Zentrum Dresden- Rossendorf, Dresden 01328, Germany Here, ϵ0 is the vacuum permittivity. Using our χeff(3) of −3 × 10−8 m2/ Mischa Bonn − Max-Planck-Institut für Polymerforschung, V and n3f = 10, we obtain an n2 close to −10 cm /GW. This number 2 6 2 Mainz 55128, Germany; orcid.org/0000-0001-6851- is valid for a very small incident field strength. We note that even with 8453 a field strength of 5 kV/cm (equivalent to an intensity of ITHz = cϵ0n| Frank H.L. Koppens − ICFO - Institut de Ciències Fotòniques, E|2/2 ≈ 10−5 GW/cm2), we do not meet the criterion above. The Barcelona Institute of Science and Technology, Therefore, we make an alternative estimate of n2 that is valid for the experimentally used field strengths. We consider a field strength of Castelldefels (Barcelona) 08860, Spain; ICREA - Institució ∼50 kV/cm, where the third-harmonic signal is strongly saturated, Catalana de Reçerca i Estudis Avancats, Barcelona 08010, indicating that the carrier temperature is so high that the THz Spain; orcid.org/0000-0001-9764-6120 conductivity is close to zero (see also refs 16 and 22). This means that Dmitry Turchinovich − Fakultät für Physik, Universität the refractive index has decreased from ∼10 to something close to 2.5 Bielefeld, Bielefeld 33615, Germany − the refractive index of graphite.40 We thus use n2 = (n − n0)/ITHz, Michael Gensch − Institute of Optical Sensor Systems, DLR, with (n − n0) = −7.5 and ITHz = 0.025 GW/cm2 (for a field strength Berlin 12489, Germany; Institut für Optik und Atomare of 50 kV/cm), obtaining n2 = −300 cm2/GW. So the nonlinear Physik, Technische Universität Berlin, Berlin 10623, refractive index for our experimental conditions is |n2| > 100 cm2/GW. Germany 1152 https://dx.doi.org/10.1021/acsnano.0c08106 ACS Nano 2021, 15, 1145−1154 ACS Nano www.acsnano.org Article Complete contact information is available at: et al. Terahertz-Field-Induced Insulator-to-Metal Transition in https://pubs.acs.org/10.1021/acsnano.0c08106 Vanadium Dioxide Metamaterial. Nature 2012, 487, 345−348. (14) Fan, K.; Hwang, H. Y.; Liu, M.; Strikwerda, A. C.; Sternbach, Notes A.; Zhang, J.; Zhao, X.; Zhang, X.; Nelson, K. A.; Averitt, R. D. Nonlinear Terahertz Metamaterials via Field-Enhanced Carrier The authors declare no competing financial interest. Dynamics in GaAs. Phys. Rev. Lett. 2013, 110, 217404. (15) Kim, S.; Jin, J.; Kim, Y.-j.; Park, I.-y.; Kim, Y.; Kim, S.-w. High- ACKNOWLEDGMENTS Harmonic Generation by Resonant Plasmon Field Enhancement. We thank David Saleta and Jake Mehew for Raman Nature 2008, 453, 757−760. characterization and analysis, Z. Wang for assistance in the (16) Hafez, H.; Kovalev, S.; Deinert, J.-C.; Mics, Z.; Green, B.; Awari, N.; Chen, M.; Germanskiy, S.; Lehnert, U.; Teichert, J.; et al. HHG experiments, and Ulf Lehnert and Jochen Teichert for Extremely Efficient Terahertz High-Harmonic Generation in Gra- technical assistance. K.-J.T. acknowledges funding from the phene by Hot Dirac Fermions. Nature 2018, 561, 507−511. European Union’s Horizon 2020 research and innovation (17) Giorgianni, F.; Chiadroni, E.; Rovere, A.; Cestelli-Guidi, M.; program under Grant Agreement No. 804349 (ERC StG Perucchi, A.; Bellaveglia, M.; Castellano, M.; Di Giovenale, D.; Di CUHL) and financial support through the MAINZ Visiting Pirro, G.; Ferrario, M.; et al. Dirac Surface States in Bi2Se3 Professorship. ICN2 was supported by the Severo Ochoa Topological Insulator. Nat. Commun. 2016, 7, 11421. program from Spanish MINECO Grant No. SEV-2017-0706. (18) Kovalev, S.; Deinert, J.-C.; Ilyakov, I.; Awari, N.; Chen, M.; Parts of this research were carried out at ELBE at the Ponomaryov, A.; Bawatna, M.; de Oliveira, T.; Eng, L. M.; Kuznetsov, Helmholtz-Zentrum Dresden-Rossendorf e.V., a member of K. A.; et al. Perturbative Regime of Terahertz High-Harmonic the Helmholtz Association. N.A., S.K., and I.I. acknowledge Generation in Topological Insulators. arXiv 2020; https://arxiv.org/ support from the European Union’s Horizon 2020 research abs/2006.03948 (accessed 2020-11-24). and innovation program under Grant Agreement No. 737038 (19) Cheng, B.; Kanda, N.; Ikeda, T. N.; Matsuda, T.; Xia, P.; Schumann, T.; Stemmer, S.; Itatani, J.; Armitage, N. P.; Matsunaga, R. (TRANSPIRE). X.J. acknowledges the support from the Max Efficient Terahertz Harmonic Generation with Coherent Acceleration Planck Graduate Center with the Johannes Gutenberg- of Electrons in the Dirac Semimetal Cd3As2. Phys. Rev. Lett. 2020, Universität Mainz (MPGC). 124, 117402. (20) Kovalev, S.; Dantas, R. M. A.; Germanskiy, S.; Deinert, J.-C.; REFERENCES Green, B.; Ilyakov, I.; Awari, N.; Chen, M.; Bawatna, M.; Ling, J.; (1) Reshef, O.; De Leon, I.; Alam, M. Z.; Boyd, R. W. Nonlinear et al. Non-Perturbative Terahertz High-Harmonic Generation in the Optical Effects in Epsilon-near-Zero Media. Nat. Rev. Mater. 2019, 4, Three-Dimensional Dirac Semimetal Cd3As2. Nat. Commun. 2020, 11, 535−551. 2451. (2) Kauranen, M.; Zayats, A. V. Nonlinear Plasmonics. Nat. (21) Soavi, G.; Wang, G.; Rostami, H.; Purdie, D. G.; De Fazio, D.; Photonics 2012, 6, 737−748. Ma, T.; Luo, B.; Wang, J.; Ott, A. K.; Yoon, D.; et al. Broadband, (3) Miller, D. A. B. Are Optical Transistors the Logical Next Step? Electrically Tunable Third-Harmonic Generation in Graphene. Nat. Nat. Photonics 2010, 4, 3−5. Nanotechnol. 2018, 13, 583−588. (4) Aouani, H.; Rahmani, M.; Navarro-Cia, M.; Maier, S. A. Third- (22) Mics, Z.; Tielrooij, K.-J.; Parvez, K.; Jensen, S. a.; Ivanov, I.; Harmonic-Upconversion Enhancement from a Single Semiconductor Feng, X.; Müllen, K.; Bonn, M.; Turchinovich, D. Thermodynamic Nanoparticle Coupled to a Plasmonic Antenna. Nat. Nanotechnol. Picture of Ultrafast Charge Transport in Graphene. Nat. Commun. 2014, 9, 290−294. 2015, 6, 7655. (5) Alam, M. Z.; De Leon, I.; Boyd, R. W. Large Optical (23) Hafez, H. A.; Kovalev, S.; Tielrooij, K.-J.; Bonn, M.; Gensch, Nonlinearity of Indium Tin Oxide in its Epsilon-near-Zero Region. M.; Turchinovich, D. Terahertz Nonlinear Optics of Graphene: From Science 2016, 352, 795−797. Saturable Absorption to High-Harmonics Generation. Adv. Opt. (6) Alam, M. Z.; Schulz, S. A.; Upham, J.; De Leon, I.; Boyd, R. W. Mater. 2020, 8, 1900771. Large Optical Nonlinearity of Nanoantennas Coupled to an Epsilon- (24) Jensen, S.; Mics, Z.; Ivanov, I.; Varol, H.; Turchinovich, D.; near-Zero Material. Nat. Photonics 2018, 12, 79−83. Koppens, F.; Bonn, M.; Tielrooij, K. Competing Ultrafast Energy (7) You, J. W.; Bongu, S. R.; Bao, Q.; Panoiu, N. C. Nonlinear Relaxation Pathways in Photoexcited Graphene. Nano Lett. 2014, 14, Optical Properties and Applications of 2D Materials: Theoretical and 5839−5845. Experimental Aspects. Nanophotonics 2018, 8, 63−97. (25) Tielrooij, K.; Song, J.; Jensen, S.; Centeno, A.; Pesquera, A.; (8) Nie, S.; Emory, S. R. Probing Single Molecules and Single Zurutuza Elorza, A.; Bonn, M.; Levitov, L.; Koppens, F. Photo- Nanoparticles by Surface-Enhanced Raman Scattering. Science 1997, excitation Cascade and Multiple Hot-Carrier Generation in Graphene. 275, 1102−1106. Nat. Phys. 2013, 9, 248−252. (9) Fan, W.; Zhang, S.; Malloy, K. J.; Brueck, S. R. J.; Panoiu, N. C.; (26) Tomadin, A.; Hornett, S.; Wang, H.; Alexeev, E.; Candini, A.; Osgood, R. M. Second Harmonic Generation from Patterned GaAs Coletti, C.; Turchinovich, D.; Kläui, M.; Bonn, M.; Koppens, F.; et al. inside a Subwavelength Metallic Hole Array. Opt. Express 2006, 14, The Ultrafast Dynamics and Conductivity of Photoexcited Graphene 9570−9575. at Different Fermi Energies. Sci. Adv. 2018, 4, No. eaar5313. (10) van Nieuwstadt, J. A. H.; Sandtke, M.; Harmsen, R. H.; (27) Goossens, S.; Navickaite, G.; Monasterio, C.; Gupta, S.; Segerink, F. B.; Prangsma, J. C.; Enoch, S.; Kuipers, L. Strong Piqueras, J. J.; Pérez, R.; Burwell, G.; Nikitskiy, I.; Lasanta, T.; Galán, Modification of the Nonlinear Optical Response of Metallic T.; et al. Broadband Image Sensor Array Based on Graphene−CMOS Subwavelength Hole Arrays. Phys. Rev. Lett. 2006, 97, 146102. Integration. Nat. Photonics 2017, 11, 366−371. (11) Genevet, P.; Tetienne, J.-P.; Gatzogiannis, E.; Blanchard, R.; (28) Green, B.; Kovalev, S.; Asgekar, V.; Geloni, G.; Lehnert, U.; Kats, M. A.; Scully, M. O.; Capasso, F. Large Enhancement of Golz, T.; Kuntzsch, M.; Bauer, C.; Hauser, J.; Voigtlaender, J.; et al. Nonlinear Optical Phenomena by Plasmonic Nanocavity Gratings. High-Field High-Repetition-Rate Sources for the Coherent THz Nano Lett. 2010, 10, 4880−4883. Control of Matter. Sci. Rep. 2016, 6, 22256. (12) Cai, W.; Vasudev, A. P.; Brongersma, M. L. Electrically (29) Ganeev, R. A.; Kulagin, I. A.; Ryasnyansky, A. I.; Tugushev, R. Controlled Nonlinear Generation of Light with Plasmonics. Science I.; Usmanov, T. Characterization of Nonlinear Optical Parameters of 2011, 333, 1720−1724. KDP, LiNbO3 and BBO Crystals. Opt. Commun. 2004, 229, 403−412. (13) Liu, M.; Hwang, H. Y.; Tao, H.; Strikwerda, A. C.; Fan, K.; (30) Vermeulen, N.; Castello-Lurbe, D.; Khoder, M.; Pasternak, I.; Keiser, G. R.; Sternbach, A. J.; West, K. G.; Kittiwatanakul, S.; Lu, J.; Krajewska, A.; Ciuk, T.; Strupinski, W.; Cheng, J.; Thienpont, H.; Van 1153 https://dx.doi.org/10.1021/acsnano.0c08106 ACS Nano 2021, 15, 1145−1154 ACS Nano www.acsnano.org Article Erps, J. Graphene’s Nonlinear-Optical Physics Revealed through Exponentially Growing Self-Phase Modulation. Nat. Commun. 2018, 9, 2675. (31) Sirtori, C.; Capasso, F.; Sivco, D. L.; Cho, A. Y. Giant, Triply Resonant, Third-Order Nonlinear Susceptibility χ(3) 3ω in Coupled Quantum Wells. Phys. Rev. Lett. 1992, 68, 1010−1013. (32) Markelz, A. G.; Gwinn, E. G. Nonlinear Response of Quantum- Confined Electrons in Nonparabolic Subbands. J. Appl. Phys. 1996, 80, 2533−2535. (33) Kono, J.; Su, M. Y.; Inoshita, T.; Noda, T.; Sherwin, M. S.; Allen, S. J.; Sakaki, H. Resonant Terahertz Optical Sideband Generation from Confined Magnetoexcitons. Phys. Rev. Lett. 1997, 79, 1758−1761. (34) Zanotto, S.; Tredicucci, A. Universal Lineshapes at the Crossover between Weak and Strong Critical Coupling in Fano- Resonant Coupled Oscillators. Sci. Rep. 2016, 6, 24592. (35) Graham, M. W.; Shi, S.-F.; Ralph, D. C.; Park, J.; McEuen, P. L. Photocurrent Measurements of Supercollision Cooling in Graphene. Nat. Phys. 2013, 9, 103−108. (36) Turchinovich, D.; Dijkhuis, J. I. Performance of Combined < 100 > − < 110 > ZnTe Crystals in an Amplified THz Time-Domain Spectrometer. Opt. Commun. 2007, 270, 96−99. (37) Yang, Y.; Lu, J.; Manjavacas, A.; Luk, T. S.; Liu, H.; Kelley, K.; Maria, J.-P.; Runnerstrom, E. L.; Sinclair, M. B.; Ghimire, S.; Brener, I. High-Harmonic Generation from an Epsilon-near-Zero. Nat. Phys. 2019, 15, 1022−1027. (38) Kovalev, S.; Green, B.; Golz, T.; Maehrlein, S.; Stojanovic, N.; Fisher, A. S.; Kampfrath, T.; Gensch, M. Probing Ultra-Fast Processes with High Dynamic Range at 4th-Generation Light Sources: Arrival Time and Intensity Binning at Unprecedented Repetition Rates. Struct. Dyn. 2017, 4, 024301. (39) Yeh, K.-L.; Hoffmann, M. C.; Hebling, J.; Nelson, K. A. Generation of 10 μJ Ultrashort Terahertz Pulses by Optical Rectification. Appl. Phys. Lett. 2007, 90, 171121. (40) Blake, P.; Hill, E. W.; Castro Neto, A. H.; Novoselov, K. S.; Jiang, D.; Yang, R.; Booth, T. J.; Geim, A. Making Graphene Visible. Appl. Phys. Lett. 2007, 91, 063124. 1154 https://dx.doi.org/10.1021/acsnano.0c08106 ACS Nano 2021, 15, 1145−1154