This document is the Accepted Manuscript version of a Published Work that appeared in final form in Nano Letters, copyright © American Chemical Society after peer review and technical editing by the publisher. To access the final edited and published work see: https://dx.doi.org/10.1021/acs.nanolett.8b04171. Fast and sensitive terahertz detection using an antenna-integrated graphene pn-junction Sebastián Castilla,† Bernat Terrés,† Marta Autore,‡ Leonardo Viti,¶ Jian Li,§ Alexey Y. Nikitin,k,⊥ Ioannis Vangelidis,# Kenji Watanabe,@ Takashi Taniguchi,@ Elefterios Lidorikis,# Miriam S. Vitiello,¶ Rainer Hillenbrand,‡,⊥ Klaas-Jan Tielrooij,∗,†,4 and Frank H.L. Koppens∗,†,∇ †ICFO - Institut de Ciències Fotòniques, The Barcelona Institute of Science and Technology, Castelldefels (Barcelona) 08860, Spain ‡CIC nanoGUNE, E-20018 Donostia-San Sebastian, Spain ¶NEST, CNR - Istituto Nanoscienze and Scuola Normale Superiore, 56127 Pisa, Italy §State Key Laboratory of Analytical Chemistry for Life Science, School of Chemistry and Chemical Engineering, Nanjing University, Nanjing 210023, China kDonostia International Physics Center (DIPC), Donostia-San Sebastin, Spain ⊥IKERBASQUE, Basque Foundation for Science, 48013 Bilbao, Spain #Department of Materials Science and Engineering, University of Ioannina, GR-45110 Ioannina, Greece @Advanced Materials Laboratory, National Institute for Material Science, 305-0044, Tsukuba, Japan 4Current address: Catalan Institute of Nanoscience and Nanotechnology (ICN2), Barcelona Institute of Science and Technology, Campus UAB, Bellaterra, Barcelona, 08193, Spain ∇ICREA - Institució Catalana de Reçerca i Estudis Avancats, 08010 Barcelona, Spain E-mail: klaas.tielrooij@icn2.cat; frank.koppens@icfo.eu 1 Abstract Although the detection of light at terahertz (THz) frequencies is important for a large range of applications, current detectors typically have several disadvantages in terms of sensitivity, speed, operating temperature, and spectral range. Here, we use graphene as photoactive ma- terial to overcome all of these limitations in one device. We introduce a novel detector for tera- hertz radiation that exploits the photo-thermoelectric effect, based on a design that employs a dual-gated, dipolar antenna with a gap of ∼100 nm. This narrow-gap antenna simultaneously creates a pn-junction in a graphene channel located above the antenna, and strongly concen- trates the incoming radiation at this pn-junction, where the photoresponse is created. We demonstrate √ that this novel detector has excellent sensitivity, with a noise-equivalent power of 80 pW/ Hz at room temperature, a response time below 30 ns (setup-limited), a high dy- namic range (linear power dependence over more than 3 orders of magnitude) and broadband operation (measured range 1.8 – 4.2 THz, antenna-limited), which fulfills a combination that is currently missing in the state of the art. Importantly, based on the agreement we obtain between experiment, analytical model, and numerical simulations, we have reached a solid understanding of how the PTE effect gives rise to a THz-induced photoresponse, which is very valuable for further detector optimization. √ Photodetectors operating at THz frequencies with an NEP on the order of ∼1 nW/ Hz. play an important role in many applications However, their response time is very long: 100 in the fields of medicine, security, quality test- and 30 ms, respectively. 11,12 Bolometric THz ing, chemical spectroscopy and more. 1–7 One of detectors, on the other hand, can√be highly the main benefits of THz radiation is its non- sensitive with an NEP of ∼0.5 pW/ Hz, while invasive nature and its capability to penetrate simultaneously showing fast operation with a most dielectric materials, which are typically response time of ∼50 ps. However, these detec- opaque at non-THz frequencies. For example, tors require cryogenic temperatures (∼4 K) and in the case of medical imaging and security suffer from a narrow dynamic range (maximum applications, THz radiation offers clear advan- detectable power ∼0.1 µW). 13 Schottky diodes, tages since it is non-ionizing due to its low although combining high speed (response time photon energy (in the meV range) in contrast in the picoseconds regime) √ and high sensitiv- with conventional X-ray radiation with much ity (NEP of 10 - 100 pW/ Hz), have a low higher photon energy (in the keV range), lead- frequency cut-off (operation only below ∼1 ing to strongly reduced health risks. Further- THz) and a small dynamic range. 14,15 Thus, more, THz detectors are expected to play an currently there are no commercially available enabling role for data communication at THz THz detectors that simultaneously meet all five bit rates. 8–10 For many of these applications, requirements. the ideal THz detector would meet the following five requirements: it should be highly sensitive Owing to its exceptional optoelectronic prop- (i.e. have a low noise-equivalent power, NEP), erties and broadband absorption spectrum operate at room temperature, give a fast pho- (from the visible down to the GHz-THz range) toresponse, have a high dynamic range (the graphene is a highly interesting photoactive range between the lowest and highest measur- material for detecting light. 16–21 During the able incident light power), and work over a past couple of years, there were several experi- broad range of THz frequencies, in particular mental demonstrations of graphene-based pho- above 1 THz. todetection in the GHz-THz range. These de- tectors were based on various operating mech- Commercially available room-temperature anisms. First of all, there were reports de- THz detectors, for example pyroelectric detec- scribing plasma wave-assisted THz detection, tors and Golay cells, are reasonably sensitive typically in the overdamped regime. 20,22–25 This 2 scheme has resulted in √ high sensitivities, with dominant photoresponse mechanism. We sup- an NEP below nW/ Hz, but typically this port this with a quantitative comparison of the mechanism is reported for radiation below 1 device response with numerical simulations and THz. Secondly, ballistic graphene rectifiers an analytical model of the PTE photoresponse. were demonstrated with excellent sensitivity, We furthermore show that, owing to its novel but only operating below 1 THz. 26 Moreover, device design, our PTE THz photodetector a graphene-antenna coupled bolometer for de- meets all five requirements of an ideal detector. tecting GHz radiation was shown with promis- In addition, it has the advantage of being based ing values of sensitivity at low temperatures. 27 on low-cost materials with scalable integration capabilities with the well established CMOS Alternatively, one can exploit the photo- electronics for low-cost imaging systems. 34 Fi- thermoelectric (PTE) effect, where absorbed nally, it is very low in power consumption, as THz light heats up the graphene electrons, 28 it is a passive device. subsequently creating an electron-heat driven photoresponse if an asymmetry is present in In the following, we will first explain how our the device. 29 Such an asymmetry could be cre- antenna-integrated pn-junction THz detector ated, for example, by using two different con- works, followed by the experimental charac- tact metals or by using two adjacent graphene terization of the detector. Subsequently, we regions of different doping, e.g. forming a junc- provide an analytical model of the PTE detec- tion. Photodetection based on the PTE ef- tor and numerical simulations of the absorp- fect in graphene was first shown for visible tion enhancement of graphene induced by the light, where interband absorption of light oc- antenna, and compare these results with the curs. 30,31 More recently, also photoresponses experiments. Finally, we will compare the THz in the THz frequency regime, where absorp- photodection performance to the state of the tion occurs through intraband processes, were art. attributed to the PTE effect. 32,33 Moreover, several of the previously mentioned studies ex- Our THz photodetector is based on a novel ploiting alternative mechanisms also attributed design (see Fig. 1a-c), which works as follows. a potentially significant fraction of the observed The detector contains a dipole antenna that THz photoresponse to the photo-thermoelectric is located ∼15 nm below a graphene channel. effect. 20–26 Clearly, some controversy exist on The antenna consists of two branches that are which photoresponse mechanism dominates for separated by a very narrow gap, with a size on graphene excited by THz light, which has ham- the order of 100 nm. This antenna gap serves pered the development of more optimized de- for focusing the incoming THz radiation into a tectors. Furthermore, the main challenge for very small spot in the graphene channel. Here, exploiting the PTE effect for THz detection is the concentrated field of the antenna leads the large mismatch between the large area of to (intraband) absorption and the subsequent the incoming radiation and the small photo- creation of hot carriers. 28 Since the creation active area of graphene, where the PTE effect of a photoresponse from hot carriers requires occurs. a gradient in the Seebeck coefficient, we use the antenna branches simultaneously as split Here, we solve this issue by introducing an gates. We apply appropriate voltages (VL and antenna-integrated THz photodetector, based VR ) to the left and right antenna branch, and on high-mobility, gate-tunable, hexagonal BN through capacitive coupling this creates a pn- (hBN)-encapsulated graphene, where the in- junction in the graphene channel, and thereby coming THz radiation is concentrated such a THz-induced photoresponse. Thus, the an- that it overlaps with the small photoactive area tenna simultaneously creates the photoactive of the graphene. Using the gate-tunability of area in the graphene channel (located around the detector, we find that the PTE effect is the the pn-junction, see Fig. 1c) and funnels inci- 3 Figure 1: a) Schematic representation (right; not to scale) of the antenna-integrated pn-junction device and a zoom of the central part of the THz PTE detector (left; to scale), consisting of an “H-shaped” graphene channel, contacted by source and drain electrodes. Underneath the graphene channel, there are two antenna branches that concentrate the incident THz light around the antenna gap region. The color map superimposed on the device shows the simulated power profile (|E/E0 |2 , where E0 is the incident electric field) at a position 5 nm below the graphene channel. The black scale bar corresponds to 1.6 µm. b) Side view of the device design, with the superimposed color map again indicating the normalized power profile as in panel a. The region where the field is strongly enhanced by the antenna overlaps with the central part of the graphene channel. c) Same as panel a, now indicating how the antenna branches serve as local gates by applying voltages VL and VR . Appropriate voltages will create a pn- junction in the central part of the graphene channel, directly above the antenna gap (which is where incident THz light is concentrated by the antenna). The color map superimposed on the device is a simulation that shows the photoresponse created by local photoexcitation, varying the position of photoexcitation (see also SI). The largest photoresponse is created when photoexcitation occurs around the junction region. The photoresponse then decreases exponentially when moving away on both sides from the junction, with the exponential decay length given by the cooling length `cool . The photoactive area (dashed rectangle) therefore has a length 2 · `cool , and a width w, which is the width of the central part of the graphene channel. d) Photocurrent image (log-scale) obtained by scanning the detector in the focal plane of a focused laser beam at 3.4 THz. We use our QCL with an average power of 84.1 µW, and a peak irradiance in the center of the focus of 1200 W/m2 . The THz light is polarized parallel to the antenna axis. The observation of the Airy pattern with multiple observable rings indicates excellent detector sensitivity. dent radiation to this photoactive area, due to cal connection between the antenna and the the very strong field enhancement of incident graphene, which means that there is no need THz radiation above the gap between the two for impedance matching to assure current flow antenna branches (see Fig. 1a-b). between antenna and graphene. The fact that we simultaneously use our antenna for focus- Compared to previous antenna-integrated, ing light and as split gate, has the advantage graphene-based THz detectors, 20–27 the ad- that there automatically is very good overlap vantage of our design is that the antenna between the region where the incoming THz gap is much smaller (100 nm vs. several mi- radiation is focused and the photoactive region crons), which means that the THz intensity of the graphene channel (see Fig. 1a-c). Fur- is greatly enhanced (∼4 orders of magnitude, thermore, we use hBN-encapsulated graphene, see Fig. 1a). Also, there is no direct electri- which leads to graphene with high mobility and 4 low residual doping. This means that the re- a pre-amplifier and/or lock-in amplifier. We sistance of the graphene channel will be low, use typical incident THz powers in the range of and we can tune the system relatively close to several microwatts to several milliwatts. the Dirac point (point of lowest carrier density), where the Seebeck coefficient is largest. Finally, We first scan our photodetector (Device A) we pattern the graphene channel into an “H- across the THz focus using motorized stages shape” with a relatively narrow (micron-sized) (in the QCL setup, see Methods). The dipolar width. The narrow width of the central part antenna, with a length of 60 µm, is smaller of the channel leads to an enhanced photore- than the THz focus (FWHM ∼200 µm), al- sponse, as the hot carriers will have a higher lowing us to spatially map out the intensity of temperature. The wider sides of the graphene the THz focus through the photocurrent IPTE . channel reduce overall device resistance by min- The results show a clear Airy pattern (see Fig. imizing contact resistance. 1d), where we are able to observe several rings of the diffracted beam pattern. This suggests We have fabricated two “H-shaped”, high- that our THz photodetector is very sensitive, mobility, hBN-encapsulated graphene devices considering that these rings contain only a very with a dipolar antenna/gating structure, (see small fraction of total incident power of the Methods and Supporting Information for de- THz beam (Pin = 84.1 µW). tails and optical pictures of the devices). We mainly show the results obtained from THz Before proceeding with quantifying the sensi- photodetector Device A and will mention some tivity, we first exploit the gate-tunability to results from THz photodetector Device B. Both identify the photocurrent generation mech- Device A and Device B have a width of the anism and determine the optimal operating central part of the graphene channel (at the point of our THz detector, by mapping out the junction) of w = 2 µm, whereas the gap sizes photoresponse as a function of gate voltages VL of the dipolar antennas are 200 and 100 nm, and VR (see Fig. 2a). These measurements were respectively. The vertical distance between an- done with the THz gas laser at 2.52 THz, and tenna and graphene is given by the thickness of using Device A. The gate voltages indepen- the bottom hBN layer, typically ∼15 nm, and dently control the carrier density (Fermi en- thus small enough to warrant sufficient overlap ergy) of the two graphene regions and therefore between the electric field profile around the an- the Seebeck coefficients S1 and S2 . The Seebeck tenna gap, and the graphene (see also Fig. 1b). coefficient of graphene has a non-monotonous The graphene mobility for both devices is on dependence on carrier density, where it first the order of 20,000 cm2 /Vs, which is a lower increases upon approaching the Dirac point bound as it is determined from two-terminal and then changes sign when crossing the Dirac measurements (see Supporting Information, point, i.e. when going from hole to electron SI). We characterize the performance of the doping or vice versa (see Fig. 2b). Since the THz photodetector devices using two different generated photocurrent IPTE ∝ (S1 − S2 ), this setups, both containing a THz laser and opti- leads to the characteristic sixfold pattern, first cal components to focus the light at our THz shown in Ref. 31 for visible light, and explained detector. One setup contains a pulsed quan- in Fig. 2b. The fact that we also observe a tum cascade laser (QCL) operating at 3.4 THz, sixfold pattern strongly suggest that our THz and two Picarin (tsupurica) lenses to focus the photoresponse is dominated by the PTE effect. light; the other setup contains a continuous To further confirm that the PTE mechanism wave THz gas laser with tunable output fre- dominates over alternative photocurrent mech- quency, and a parabolic mirror to focus the anisms, such as bolometric and photogating light (see Methods for details). The THz light effects, we measured the photocurrent as a func- is usually modulated by an optical chopper and tion of bias voltage applied between the source the generated photocurrent is measured using and drain contacts. The drain current increases 5 Figure 2: a) Photoresponse as a function of voltages applied to the two antenna branches/gates, with radiation at 2.52 THz. We use our THz gas laser with an average incident power of Pin = 5 mW, and a peak irradiance in the center of the focus of 1.8·104 W/m2 . The THz light is polarized parallel to the antenna axis. The sixfold pattern, as in Refs., 29,31 indicates that the photoresponse is generated through the PTE effect. The photoresponse is the photocurrent normalized by the power in a diffraction-limited spot, i.e. the responsivity R. The maximum responsivity occurs in the pn- and np-regions. b) Line cut at the location of the dashed line in panel a, showing a double sign change (red dots and line; left vertical axis) as a function of carrier density (controlled through gate voltage VL ). The blue line represents the Seebeck coefficient (calculated from the experimentally obtained graphene mobility, see SI; right vertical axis). The double sign change occurs due to the non-monotonous dependence of the Seebeck coefficient on carrier density: for a constant Seebeck coefficient in one region (dashed horizontal line), the Seebeck coefficient of the other region is first higher, then lower and then again higher, giving rise to two sign changes, as indeed observed experimentally. c) The extracted Johnson noise current, calculated from the resistance that was measured simultaneously with the result in panel a. d) The noise-equivalent power (NEP), extracted from the results in panels a and c, normalized to a diffraction-limited spot. The white cross indicates the gate configuration that corresponds to the lowest NEP: the left (right) gate at 0.20 V (-0.72 V), corresponding to an electron density of 7.5·1011 cm−2 , EF = +100 meV (hole density of 3.6·1011 cm−2 , EF = −70 meV). linearly with applied bias voltage, whereas the photocurrent remains constant (see SI), in con- We now proceed with quantifying the sen- tradiction with what is expected for the bolo- sitivity of our THz photodetector. First, we metric and photogating effects. Thus, these identify the largest responsivity R = IPTE /Pdiff results show that the photo-thermoelectric ef- at the optimal gate configuration for Device A. fect is responsible for the observed THz pho- We note that the responsivity that we use is toresponse. We find the largest photoresponse the responsivity normalized by the power in a in the pn-junction and np-junction regimes, as diffraction-limited spot with NA = 1. In some expected, relatively close to the Dirac point. of the literature, a responsivity is provided that 6 is normalized by the amount of power that is a DevA A Dev. actually incident in the experiment Pin . How- -7 1.83 2.52 100 10 3.11 ever in the case where the wavelength is larger 4.25 IPTE [nA] Dev. DevB B 2.52 THz than the photodetector device this number does 10 10 -8 not only depend on the device performance, but also on how well the THz light is focused. 1 10 -9 Alternatively, one can use the responsivity nor- malized to the power that is absorbed in the 0.1 10 -7 1 10 -6 10 -5 10 100 -4 10 actual graphene channel or in the photoac- Pdiff [µW] tive area, which would yield an artificially high b 15 1200 enhancement G number, 32 as it is impossible to focus the THz 6 6 Antenna  [mA/W] light in such a small area. Here we choose the 10 800 responsivity normalized by the incident power 4 4 in a diffraction-limited spot, because this is ar- 5 400 2 2 guably the most technologically relevant num- ber (as this represents what will be reached 0 0 1 1 2 2 3 3 4 4 5 5 when combining the detector with an optimized Frequency [THz] focusing system, such as using a silicon hemi- c 0.2 spherical lens 23–25 ), and it is the convention 0.0 -0.1 laser that is most commonly used in the literature 0.0 -0.2 off on VPTE [V] -0.3 on THz photodetection (see also Table 1). We -0.2 -0.4 00 Time10[µs] 15 5 20 20 calculate the power in a diffraction-limited spot using Pdiff = Pin · Adiff /Afocus , where Pin is the -0.4 -0.4 measured total incident THz power, Adiff is the -0.6 calculated area of a diffraction-limited spot and 0 50 100 150 200 Afocus is the measured area of the focused THz Time [ns] beam. Typically, we have Adiff /Afocus ≈ 1/60 (corresponding to a NA of ∼0.13 for our fo- Figure 3: a) Photocurrent as a function of THz laser cusing system based on a parabolic mirror, see power (in a diffraction-limited spot) Pdiff in log-log scale. All round data points correspond to Device A, Methods for details). From Fig. 2a we extract whereas the light blue squares correspond to Device B. a maximum responsivity of R = 14 mA/W (32 The lines through the experimental data points are fits V/W). For Device B, we find a maximum re- according to IPTE ∝ Pdiff γ . The obtained exponent is sponsivity value of R = 25 mA/W (105 V/W) close to 1 for all data sets. For Device A, these data (see Methods and SI). In both cases the THz correspond to a gate configuration of VL = 0 V and VR = -0.67 V, corresponding to an electron (hole) den- light was at 2.52 THz (corresponding to a wave- sity of 4.2·1011 (-2.7·1011 ) cm−2 . b) Responsivity R length of 118.96 µm, 84 cm−1 ). as a function of THz wavelength (red dots, left vertical axis), with the same (sub-optimal) gate configuration as in panel a. The black line (right vertical axis) shows the result for the antenna-induced absorption enhancement in the graphene channel. The red dashed line illustrates the trend of the experimental points. c) Results of the pulsed laser experiment, where the photocurrent was amplified by a fast current pre-amplifer (Femto) and the data were acquired with a fast oscilloscope. The inset shows how the photovoltage VPTE follows the switching of the pulsed laser. The red and blue (plotted with an offset) open dots show the obtained photovoltage in a small time window marked in the inset, with the black line giving the result of exponential fits with timescales of 40 (24) ns for the red (blue) curve. 7 Using the extracted responsivity, we now de- expect a sub-linear dependence of photocurrent termine the sensitivity of the detector. For this, on power, with an exponent that tends to γ = we note that our detector operates without bias, 0.5. which means that it is limited by theqJohnson or thermal noise, given by: Inoise = 4kBRT ∆f . Regarding the range of frequencies where our Here, kB is the Boltzmann constant, T is the detectors operate, we note that this is only temperature of operation, ∆f is the spectral limited by the antenna structure. The rea- bandwidth and R the resistance. We show the son behind the spectrally ultra-broad photore- Johnson noise current of Device A as a function sponse of graphene is the efficient heating of of the two gates voltages in Fig. 2c, calculated the electrons, which occurs irrespective of the from the measured resistance R. The resis- wavelength of the incident light, i.e. whether tance was measured simultaneously with the intraband or interband light absorption oc- photoresponse in Fig. 2a, which we used to curs. 28,35 We characterize the spectral response determine the responsivity R. As expected, of our detector by measuring the responsivity, we see that a higher noise current occurs for while varying the frequency from 1.83 to 4.25 lower resistance values (away from the Dirac THz (see Fig. 3b). We observe a trend where point), whereas we obtain lower noise values the responsivity peaks around 3 THz. This closer to the Dirac point, where the graphene corresponds reasonably well with the antenna resistance is higher. The photodetection sen- being optimized for a frequency of 2 THz using sitivity is given by NEP= Inoise /R, which we full wave simulations (see Methods and Fig. show as a function of the two gate voltages in 3b). The discrepancy likely comes from the Fig. 2d. The lowest values for the NEP occur fact that a simplified structure was simulated, at the pn- and np-regions close to the Dirac which didn’t contain all the metallic parts that point, where the responsivity is highest and the actual device has. Importantly, these re- the noise is lowest. The sults confirm that the spectral range where the √ lowest value of the THz detector operates is currently limited only NEP map √ was 200 pW/ Hz for Device A and 80 pW/ Hz for Device B (see Methods and SI). by the antenna. Thus, using more spectrally broad antennas or a combination of antennas An important characteristic of an ideal THz one could extend the spectral range of our pho- detector is having a large range of powers over todetector, covering the spectrum all the way which the response is linear, i.e. a large dynamic from the ultraviolet, through the visible and range. Thus, we measure the photocurrent vs. infrared to the terahertz. Pdiff for four different THz frequencies as shown in Fig. 3a. We vary the THz power over more In the following, we discuss the speed of our than 3 orders of magnitude (using the THz PTE THz photodetector. We analyze the speed gas laser setup), and fit the data with a simple in Fig. 3c, where we rapidly switch the THz ra- power law IPTE ∝ Pdiff γ . We obtain for Device A diation on (white) and off (yellow) using our at 2.5 THz a power of γ = 1.1 ± 0.2 and for De- pulsed THz QCL (see Methods). We observe vice B γ = 1.1±0.15 (95% confidence intervals). that the photoresponse of our detector VPTE This shows that the photoresponse depends lin- closely follows the laser switching behavior. early on the THz power over a range of more We quantify the detector speed by fitting the than three orders of magnitude. The reason for VPTE rise and fall dynamics (see also SI) with the linear photoresponse as a function of power simple exponential equations, and obtain an is the fact that the photodetector operates in exponential (1/e) response time of 32±11 ns, the weak heating regime, where ∆T  Tambient , corresponding to a bandwidth of ∼5 MHz. In i.e. the change in temperature of the electronic this measurement, the speed is most likely lim- system is smaller than the ambient tempera- ited by the measurement electronics, namely ture Tambient . When ∆T approaches Tambient , we by the 3.5 MHz bandwidth of the current pre- amplifier, rather than by the PTE THz detector 8 itself. Indeed, the intrinsic speed of the detec- ∆T can be large in graphene (up to several tor is expected to be significantly higher, since thousand K), because of efficient heating of ultimately the intrinsic response time is limited the electrons after absorbing THz light, due by the RC-time of the detector. Performing to strong electron-electron interactions, and this calculation (see Methods), we obtain a rise because the hot carriers are relatively weakly time of 56 ps, corresponding to a detection coupled to the crystal lattice. 28 rate of 6 GHz for our device with a mobility of 20,000 cm2 /Vs. For a device with a mobility Our photodetector design maximizes the PTE of 100,000 cm2 /Vs (see SI), we find 9 ps, cor- THz photoresponse, particularly by maximiz- responding to 40 GHz. Thus, extremely fast ing ∆T and minimizing R. From a simple heat THz photodetection with switching times in equation, the temperature increase ∆T (aver- the picosecond range should be possible. aged over all charge carriers in the photoactive area) is given by We now discuss more of the underlying physics of the PTE detector, using a simple an- Pabs ∆T ≈ (2) alytical model that provides the rationale be- Aactive Γcool hind our detector design. Owing to the differ- where Pabs is the amount of THz power that ence in Seebeck coefficients at the pn-junction, is absorbed in the active area of the graphene a local photo-thermoelectric voltage is created, channel and Γcool is the heat conductivity that which leads to the flow of a photocurrent be- describes the coupling of the heated electron tween the source and drain contacts that are systems to its environment. The photoactive connected to the graphene channel. The PTE area is given by Aactive = 2`cool · w, where `cool photocurrent is then given by 29 is the hot-carrier cooling length, which can be seen as the length scale over which hot car- (S1 − S2 ) ∆T riers can move before cooling down (typically IP T E = , (1) R 0.5 – 1 µm at room temperature, 29,31,35,37,38 see where S1 and S2 are the Seebeck coefficients also Fig. 1c). In the case of hBN-encapsulated (also called thermopower) of the two regions of graphene, Γcool is the out-of-plane, interfacial the graphene channel that are independently heat conductivity where hot graphene carriers controlled by the gates/antenna branches, ∆T couple to hyperbolic hBN phonons. 37 We opti- is the temperature increase of the electronic mize ∆T by maximizing Pabs and minimizing system induced by THz radiation, and R is the Aactive . We maximize the amount of absorbed total electrical resistance, accounting for the THz power Pabs by using a dipole antenna with graphene and contact resistances. It is worth a narrow gap, which focuses the incoming THz mentioning that this current is generated under radiation down to the small (compared to the zero applied source-drain bias voltage, resulting THz radiation wavelength) graphene photoac- in very low detector noise (Johnson noise) and tive region. We further maximize ∆T by using extremely low power consumption. a narrow channel width w of 2 µm. Basically, the smaller the area where the incident power Graphene is an ideal material to exploit the is absorbed, the smaller the amount of elec- PTE effect for THz detection, because the trons that will share the heat, and therefore term (S1 − S2 ) ∆T can be large and R is typi- the larger the increase in temperature of the cally small, in particular for high-quality, hBN- electronic system, ∆T . encapsulated graphene. The Seebeck coefficient of graphene is intrinsically quite large, on the In order to further increase the responsivity, order of 100 µV/K 36 and S1 and S2 are inde- we reduce the overall resistance R of the device, pendently tunable through the gates, meaning and optimize the shape of the graphene chan- that (S1 − S2 ) can be maximized. Furthermore, nel. We achieve low R by using high-quality hBN-encapsulated graphene (see Fig. 1b). 39,40 9 This method enables mobility values as high radiance as in the experiment (see Methods). as 100,000 cm2 /Vs at room temperature, and These simulations give an absorbed power (in low levels of intrinsic doping (see SI). We fur- the entire graphene channel) of Pabs,sim = 7 thermore pattern the graphene channel in an µW (at 2.5 THz). This number is close to the “H-shape” (see Fig. 1a), in order to reduce number we obtained experimentally, adding the overall device resistance R. This shape is credibility to our assignment of the PTE as the crucial because it has a small width w in the dominant photoresponse mechanism and to the central part of the channel – ensuring small validity of our analytical model. We ascribe the Aactive and thus large ∆T – while having a lower experimental value (by a factor ∼3) to larger width towards the contacts – minimizing non-optimal performance of the actual antenna the graphene sheet resistance. Furthermore, in the photodetector device, most likely due the large interface with the source and drain to the presence of metallic regions around the contacts minimizes the contribution of contact antenna (see SI). Furthermore, the simulations resistance to the overall resistance R. We as- consider the absorption in the entire graphene sess the validity of the analytical model of Eqs. sheet, rather than only in the photoactive area 1-2, by comparing the results with numerical of the graphene channel. Notably, we point out simulations of the PTE photocurrent gener- that without the antenna, the amount of inci- ated in different graphene geometries (see SI). dent THz light from a diffraction-limited spot We find agreement between the analytical and that would be absorbed in the photoactive area numerical results, showing the validity of our of the graphene channel would be more than analytical approach. Importantly, the analyti- three orders of magnitude lower, highlighting cal model gives us insights into the physics that the importance of the antenna-integration. We determines the detector response, thus allowing illustrate this in Fig. 3b, where we show the for optimization strategies. The advantage of antenna-induced absorption enhancement, de- Pabs,sim,w/ antenna the numerical simulations is that they are also fined as G = Pabs,sim,w/o . These simulations antenna valid for non-rectangular graphene shapes. Our show that the antenna enhances the graphene detector design is the result of these analytical absorption by more than three orders of mag- and numerical simulations. nitude. Based on the analytical model for the PTE Finally, we compare the performance of our response, we now examine our experimental photodetector with respect to other graphene- results quantitatively. We have measured a based THz photodetectors in the literature. photocurrent of IPTE = 1.14 µA (Device A, 2.5 We first compare with other detectors where THz) for an incident power of Pin = 5 mW (fo- the photodetection mechanism was explicitly cused to a spot size Afocus , see Methods). Using ascribed to the PTE effect (see Table 1). Since Eq. 1 with (S1 − S2 ) = 160 µV/K (estimated not every report used the same power normal- from Ref. 41 ) and R = 2.3 kΩ (measured), we ization procedure for the responsivity and NEP, find an experimental temperature increase of we mention explicitly the normalization proce- ∼16 K (confirming the weak heating regime). dure that was used. We note that our THz Then using Eq. 2 with interfacial heat con- detector is 2–4 orders of magnitude more sen- ductivity Γcool = 7·104 W/m2 K (determined in sitive than any other THz PTE photodetector Ref. 37 ), cooling length `cool = 510 nm (from the (if the same normalization procedure would be mobility and interfacial heat conductivity, see applied). We attribute this to our novel design SI), and channel width w = 2 µm (measured), with the antenna/gating structure, the opti- we find the absorbed power (in the active area mal graphene channel geometry and the use of the graphene channel) to be Pabs = 2.3 of high-mobility hBN-encapsulated graphene. µW. We compare this value with the absorbed Furthermore, the sensitivity that we obtain is power we find from numerical simulations of the very similar to, or better than, the most sen- antenna-graphene structure, using the same ir- 10 Table 1: Comparison of graphene-based THz photodetectors. Reference Mechanism NEP Normalization Speed Freq. range √ (pW/ Hz) area Anorm a (ns) (THz) This work PTE 80 λ2 /π <30 1.8 – 4.25 32 PTE 1100 ∼ λ2 /3350 0.11 2.5 33 PTE 350 ∼ λ2 /450 9000 0.08 – 0.3 24 Rectification 130 none - 0.4 22 Plasma waves 2000 λ2 /4 - 0.29 – 0.38 25 PTE/plasma waves 600 noneb - 0.13 – 0.45 26 Ballistic rectification 34 λ2 /4π - 0.07-0.69 a The normalization area Anorm refers to the area to which the incident power was normalized: Pnorm = Pin · Anorm /Afocus . In our work, for example, we use the power in a diffraction-limited spot Pdiff = Pin · Adiff /Afocus , i.e. we use Anorm = Adiff = λ2 /π. b Whereas the incident power was not normalized to any area, a correction of the incident power was applied to account for losses occurring in the focusing system. sitive graphene-based THz detectors reported active area, thus achieving a large THz-induced in the literature. 24,26 The operation of those change in temperature, and a low overall device detectors, however, has only been shown for resistance. frequencies below 1 THz and no response times have been measured. Additionally, it’s impor- Given the qualitative and quantitative under- tant to point out that the Drude optical con- standing we have developed of the performance ductivity and therefore absorption in graphene of our detector, we identify strategies for further is higher in the GHz range than in THz, 42,43 improving its performance. Most importantly, hence a direct comparison with detectors oper- by optimizing the antenna, a higher absorp- ating below 1 THz is not straightforward since tion and therefore lower NEP can be achieved. we do not normalize the detector responsivity Additionally, by using a more broadband an- by the graphene absorption. tenna, the detector will be sensitive for a larger range of THz frequencies. The sensitivity can In conclusion, we have demonstrated a novel be further enhanced by having a lower ther- THz photodetector, which is dominated by mal conductivity Γcool . This could be achieved the photo-thermoelectric effect. It operates at by exploring alternative encapsulation materi- room temperature, is highly sensitive and very als, rather than hBN, e.g. a transition metal fast, has a wide dynamic range and operates dichalcogenide (TMD) material, and by oper- over a broad range of THz frequencies. We have ating at a lower temperature. We estimate that optimized the PTE THz detector by using a it will √ be possible to reach an NEP in the low split-gate/antenna structure with narrow gap, pW/ Hz-regime. We expect that the unique which funnels the incident THz light exactly combination of high sensitivity and fast opera- at the small photoactive area of the detector tion means that these THz PTE detectors will leading to strongly enhanced THz absorption play an important role in a large spectrum of in graphene. This structure simultaneously applications. allows for tuning the detector to the optimal gating configuration, where a pn-junction is created in the graphene channel. Furthermore, we have used an “H-shaped”, high-quality, Supporting Information hBN-encapsulated graphene channel with a Raman, electrical characterization and optical narrow width, in order to have a small photo- images of the devices, fall time measurements of 11 device A, power dependence and scanning pho- larized (only 2% of residual unpolarized light) tocurrent measurements of device B and ther- and mounted the detector with the antenna axis moelectric simulations. parallel to the polarization. The THz light was focused using a gold parabolic mirror with focal distance 5 cm. The THz power was measured Methods using a pyroelectric THz detector from Gentec- EO placed at the sample position. The sec- Sample fabrication ond setup (used for Fig. 1d and 3c) contained a First we patterned the antenna/gate structure pulsed QCL at 3.4 THz, with an expected rise on transparent SiO2 (Infrasil) substrate using time <1 ns. electron beam lithography followed by evapo- ration of titanium (2 nm) / gold (30 nm). The Responsivity calculation antenna gap was 200 nm (100 nm) for device A For the responsivity R = IPTE /Pdiff , we ex- (B). We then released an hBN/graphene/hBN tracted the PTE photocurrent IPTE from the stack onto the antenna/gate structure. The output signal √ of the lock-in amplifier VLIA us- stack elements (top and bottom hBN and ing IPTE = 2π4G 2 VLIA , 20,22 where G is the gain graphene) were mechanically cleaved and exfo- factor in V/A (given by the lock-in amplifier). liated onto freshly cleaned Si/SiO2 substrates. The power in a diffraction-limited beam is given The full stack was prepared by the Van der by Pdiff = Pin · Adiff /Afocus , where the ratio Waals assembly technique 39,40 and released w2 onto the antenna/gate structure. This was Adiff /Afocus = w0,x0,diff w0,y . In order to obtain w0,x followed by patterning source and drain elec- and w0,y we use our observation that the pho- trodes, using electron beam lithography with a toresponse is linear in laser power and mea- PMMA 950 K resist film and exposing it to a sured the photocurrent while scanning the de- plasma of CHF3 /O2 gases for partially etching vice in the x− and y−direction. The photocur- the stack. Consequently, we evaporate side con- rent is then described by Gaussian distributions 2 2 2 2 tacts of chromium (5 nm) / gold (60 nm) and ∝ e−2x /w0,x and ∝ e−2y /w0,y , where w0,x and lift off in acetone as described in Ref. 39 . Finally, w0,y are the respectively obtained spot sizes (re- the encapsulated graphene was etched into the lated to the standard deviation via p σ = w0 /2 “H” shape using a plasma of CHF3 /O2 gases. and to the FWHM via FWHM = 2 ln(2)w0 ). From gate-dependent measurements (varying For Device A, we obtained w0,x = 263.3 µm and VL and VR simultaneously) on Device A, we ex- w0,y = 331.2 µm. For the diffraction-limited tract a mobility >20,000 cm2 /Vs and a contact spot, we took w0,diff = πλ , with λ the THz laser resistance of 126 Ω (3.8 kΩ · µm) (see SI). wavelength. The diffraction-limited area is thus 2 taken as Adiff = πw0,diff = λ2 /π. Measurements Antenna simulations In one setup (used for all Figs. except for Figs. 1d and 3c), we used a continuous wave THz The full wave simulations were performed in beam from a gas laser (FIRL 100 from Edin- Comsol. The frequency dependent permittiv- burgh Instruments) providing a maximum out- ity of hBN was taken from Ref. 44 . The op- put power in the range of a few tens of milli- tical conductivity of graphene was calculated watts, and frequencies of 1.83, 2.52, 3.11 and employing the local random phase approxima- 4.25 THz. The device position was scanned us- tion at T = 300 K with a scattering time of ing a motorized xyz-stage. The THz laser was 100 fs. In the simulations, for simplicity, the modulated at 266 Hz using an optical chop- Fermi energy of the graphene sheet was spa- per and the photocurrent was measured us- tially constant. A plane wave source was used ing a lock-in amplifier (Stanford). We verified for illumination, where the incident power was that the output THz beam was strongly po- normalized to give the same irradiance as the 12 experiment. knowledges funding from the Spanish Ministry of Economy, Industry and Competitiveness, na- Speed calculations tional project MAT2017-88358-C3-3-R. The electrically-limited operation frequency of the detector is related to the RC-time con- stant τ = RC, with R total graphene resis- tance (including contact resistance), and C to- tal graphene capacitance. The operating speed is then given by the rate f = (2πτ )−1 . The rise time τrise is the measure of the photodetector re- sponse speed to a stepped light input signal. It is the time required for the photodetector to in- crease its output signal from 10% to 90% of the final steady-state output level. The rise time is calculated as τrise = τ · ln(9) = (2πf )−1 · ln(9) = 0.35/f . Acknowledgments The authors thank Mark Lundeberg and Ia- copo Torre for fruitful discussions. F.H.L.K. acknowledges financial support from the Span- ish Ministry of Economy and Competitive- ness, through the Severo Ochoa Programme for Centres of Excellence in R&D (SEV-2015- 0522), support by Fundacio Cellex Barcelona, Generalitat de Catalunya through the CERCA program, and the Agency for Management of University and Research Grants (AGAUR) 2017 SGR 1656. Furthermore, the research leading to these results has received funding from the European Union Seventh Framework Programme under grant agreement no.785219 Graphene Flagship (Core2). ICN2 is supported by the Severo Ochoa program from Spanish MINECO (Grant No. SEV-2017-0706). K.-J.T. acknowledges support from a Mineco Young Investigator Grant (FIS2014-59639-JIN). S.C. acknowledges funding from the Barcelona In- stitute of Science and Technology (BIST), the Secretaria dUniversitats i Recerca del Depar- tament dEconomia i Coneixement de la Gen- eralitat de Catalunya and the European So- cial Fund FEDER. M.S.V. acknowledges fi- nancial support from the ERC Project 681379 (SPRINT) and partial support from the second half of the Balzan Prize 2016 in applied photon- ics delivered to Federico Capasso. A.Y.N. ac- 13 Supporting Information 1) Thermoelectric simulations To calculate the PTE response of our device we solved numerically the linearized thermoelectric equations: 45 jQ (r) = −σ(r)∇V (r) − σ(r)S(r)∇T (r), (3) jE (r) = −Π(r)σ(r)∇V (r) − [κ(r) + Π(r)σ(r)S(r)]∇T (r), (4) where jQ (r), jE (r) are the electric and energy current density respectively, V (r) is the voltage, T (r) is the local temperature, σ(r) is the electrical conductivity, S(r) is the Seebeck coefficient, Π(r) is the Peltier coefficient, and κ(r) is thermal conductivity. These two equations are coupled to the continuity equations for the charge and energy currents: ∇ · jQ (r) = 0, (5) ∇ · jE (r) = −g(r)[T (r) − T0 ] + JE (r). (6) Here g(r) parametrizes the thermal conduction to the substrate, 37 T0 is the substrate temperature, and JE (r) is a local heat source, that in our case is due to light absorption. We solved these equa- tions inside a 2D simulation box [0, Ltot ] × [0, Wtot ] illustrated in Fig. S5 that includes the graphene sheet and the gold contacts. We artificially included a small rectangular region of dimension δ × W between the graphene and the gold contacts to account for finite contact resistance by choosing the value of its conductivity σ = δ/Rcont with Rcont being the value of the gold-graphene contact resis- tance. We assumed that the material parameters σ(r), S(r), Π(r), κ(r), and g(r) are piecewise constant in the regions described in Fig. S5. Since, because of Onsager relations, Π(r) = T0 S(r), and the thermal conductivity is related to the electrical conductivity by the Wiedemann-Franz law κ(r) = L0 σ(r)T0 , with L0 ≡ π 2 kB2 /(3e2 ) = 2.44 · 10−8 WΩK−2 , we have only three independent parameters (σ(r), S(r), g(r)) that are listed in table S1. Eqs. (3,4,5,6) need to be supplemented with boundary conditions (BCs). We used Dirichlet BCs at x = 0, Ltot . These read: V (x = 0, y) = Vsd , (7) V (x = Ltot , y) = 0, (8) T (x = 0, y) = T (x = Ltot , y) = T0 . (9) We assumed instead homogeneous Neumann condition on the remaining part of the boundary n̂ · jQ (r) = 0, (10) n̂ · jE (r) = 0, (11) n̂ being the outward normal unit vector. We solved numerically Eqs. (3,4,5,6), with BCs (7,11) using the Finite Volume Method (FVM) on a regular square grid of 400x240 cells. Once the solution is found, the current flowing in the device can be calculated as Z Wtot Z Wtot Vsd I≡ Jx (0, y)dy = Jx (Ltot , y)dy = + IPTE , (12) 0 0 R where the last equality holds because of linearity, R is the resistance, and 14 Z IPTE = drRint (r)JE (r). (13) To find R we simply calculate the current by setting JE (r) = 0. The responsivity Rint (r0 ) at a given position r0 is instead calculated by setting Vsd = 0 and JE = δ(r − r0 ), and calculating the corresponding current. This step is repeated for every r0 in the simulation grid to obtain the responsivity maps in Fig. S5. Table S1: Material parameters Material σ[S] S[µV /K] Γcool [W/m2 K] Graphene n (p) doped 1.3·10−3 a 80 (-80) 41 7·104 37 Gold 4b 0c 300·104 46,47 Interface contact 2.5·10−6 d 0c 100·104 a When considering EF = 50 (-50) meV, T = 300 K and 200 fs relaxation time for both doping type n (p). b For 100 nm gold. c We neglect the contribution of the Seebeck coefficient of the metal. d For a δ = 25 nm. 15 2) Supplementary figures a b 60 μm 4 μm c d 300 μm 4 μm Figure S1: Optical microscope images of the graphene-based THz photodetectors. a) 10x magnifica- tion image of Device A, showing the metallic rods contacting the antenna branches and source-drain electrodes. b) 100x magnification image of Device A, displaying the area of the photodetector con- taining the antenna/gate structure, source-drain electrodes and “H-shaped”, h-BN encapsulated graphene. c) 5x magnification image of Device B and d) 100x magnification image of Device B. 16 Figure S2: a) Measured resistance map as a function of the right gate (left axis) and left gate (bottom axis) for Device A. b) Resistance as a function of graphene charge carrier density (linecut in dashed line shown in a)). We fit the resistance R curve using the model from Ref. 48 , where R = Rc + (L/w)(µen)−1 , Rc is contact resistance, L is channel length, pw is channel width, µ is mobility, e is elementary charge, and the carrier density is given by n = (n∗ )2 + (β(Vg − VDirac ))2 . Here n∗ is the residual doping concentration, Vg is gate voltage, VDirac is the gate voltage that corresponds to the Dirac point, and β = 0 hBN /dhBN e, where 0 is vacuum permittivity, hBN the dielectric constant of hBN and dhBN the thickness of the bottom hBN. We obtained 22,000 and 19,000 cm2 /Vs for electron (red) and hole (blue) mobility, respectively. We extracted a contact resistance of Rc = 126 Ω (3.8 kΩ · µm), and a residual doping concentration of n∗ = 1.6·1011 cm−2 . c) Carrier mobility as a function of graphene charge carrier density extracted from the conductance measurement, 39 where blue (red) corresponds to hole (electron) mobility. Figure S3: Raman spectroscopy measurements of Device A. a) Histogram plot of the full-width- half-maximum of the graphene 2D peak (Γ2D ) across the ∼15x15 µm2 region marked in red on the optical picture (see inset). The white scale bar corresponds to 30 µm. The inset shows a typical Raman spectrum of the monolayer graphene. The mean Γ2D is ∼17.5 cm−1 , characteristic of high quality single layer graphene. 49 b) 2D-peak versus G-peak frequency (ω2D and ωG , respectively) extracted from the same Raman map as in panel a. The colorbar represents the Γ2D of the recorded spectrum. The data shows low doping, as confirmed in transport measurements (see Fig. S2), and moderate levels of strain. 49 17 Figure S4: Characterization of the Hall bar reference device fabricated following the same procedure as the photodetector (see Methods). a) Optical image of the initial hBN/graphene/hBN stack. b) Optical image shown in panel a with enhanced contrast. The graphene flake is optically visible (compare the graphene shape with the Raman map shown in the inset of panel c). c) Histogram plot of the Raman scan (see inset) showing the full-width-half-maximum of the graphene 2D Raman peak (Γ2D ). The mean Γ2D is ∼ 18.1 cm−1 , characteristic of high quality single layer graphene. d) Scattering plot of the 2D-peak versus the G-peak frequency (ω2D and ωG , respectively), where the colorbar represents the Γ2D of the recorded spectrum (inset of panel c). The Raman data shows comparable behavior as the photodetector (see Fig. S3b) with low Γ2D , low doping and moderate levels of strain. 49 e) Resistance as a function of the charge carrier density of the measured Hall bar device (see inset). We fit the resistance R curve using the model from Ref. 48 , and obtained 103,000 and 96,000 cm2 /Vs for the electron (red) and hole (blue) mobilities, respectively. f ) Carrier mobility as a function of the charge carrier density extracted from the conductance measurement 39 , where blue (red) corresponds to hole (electron) mobility. These results show that a mobility of 100,000 cm2 /Vs for these kind of photodetector devices is realistic, and therefore a detection speed in the 10 ps-range can be achieved. 18 Graphene (p doped) Graphene (n doped) Metal Interface contact a b c w = 6 µm 1 µm Wtot 1 µm 6 µm Ltot Lsd = 3.6 µm 3.6 µm 3.6 µm d 0.12 e 0.12 f 0.45 Ⴆint,sim [A/W] Ⴆint,sim [A/W] Ⴆint,sim [A/W] 0.06 0.06 0.15 0.0 0.0 0.0 -0.15 -0.06 -0.06 -0.12 -0.12 -0.45 R = 3.8 kΩ R = 22.2 kΩ R = 5.3 kΩ Ⴆint,sim_max = 0.126 A/W Ⴆint,sim_max = 0.130 A/W Ⴆint,sim_max = 0.542 A/W Ⴆint,sim_tot = 0.040 A/W Ⴆint,sim_tot = 0.043 A/W Ⴆint,sim_tot = 0.48 A/W    NEP = 52.2 pW/ Hz NEP = 20.1 pW/ Hz NEP = 3.7 pW/ Hz Figure S5: Thermoelectric simulations results, 50 as in Fig. 1c, which represent the resulting pho- toresponse after local photoexcitation, expressed in photocurrent normalized by absorbed power, i.e. internal responsivity (Rint,sim as defined in Eq. S11, where JE (r) represents the locally absorbed power). These simulations were based on similar characteristics of the measured samples, with a graphene scattering time of 200 fs (mobility 20,000 cm2 /Vs), contact resistance of 10 kΩ·µm, 3.6 µm of source-drain distance (Lsd ), Γcool = 7·104 W/m2 K, EF = 50 meV and -50 meV for n- and p-doped graphene regions respectively. Panels a-c) show the simulated device geometries and panels d-f ) the internal responsivity (Rint,sim ) calculated at each position across the device. The R, Rint,sim max , Rint,sim tot and NEP below panels d-f ) indicate total device resistance (including graphene sheet and contact resistance), the maximum internal responsivity generated at the pn junction, total internal responsivity (which takes into account also the opposite sign responsivity generated at the metal-graphene interface, i.e. Rint,sim tot = Rint,sim max − | Rint,sim min |) and noise-equivalent power, respectively. The simulations show that decreasing the width of the graphene channel (going from the design of panel a to the design of panel b) does not lead to any change in responsivity, but a decrease in noise-equivalent power. The responsivity is constant because of the trade-off between the increased ∆T due to the smaller active area and the increased device resistance due to the smaller width. The NEP is reduced, because the higher device resistance leads to reduced thermal noise. Then, by using the design with an “H-shaped” graphene channel in panel c, we both increase responsivity and decrease NEP. The main reason is that the active area is significantly reduced, without increasing too much the device resistance. Furthermore, the maximum responsivity and total internal responsivity are now very similar, because the PTE response at the pn-junction dom- inates. Therefore, this is the preferred design. We compare the simulation results for the internal responsivity Rint,sim with the analytically obtained values (Rint,analyt ) using the equations in the main text and the device resistances obtained from these simulations, and obtain Rint,analyt = 0.12, 0.13 and 0.55 A/W for the designs from panels a), b) and c) respectively. These analytical values are in excellent agreement with the Rint,sim max values obtained from simulations. 19 0.0 -0.2 VPTE [V] -0.4 -0.4 -0.6 00 50 50 100 150 100 150 200 200 Time [ns] Figure S6: Results of the pulsed laser experiment, showing the fall time corresponding to the two regions (red and blue, the latter one plotted with an offset) of the inset of Fig. 3c of the main text. The open dots show the acquired photovoltage and the black lines represent the fit results with the exponential fits. We obtain an exponential fall time of 47 (28) ns for the red (blue) curve. As in the case of the rise time determination in the main text, this timescale is limited by the current amplifier with a bandwidth of 3.5 MHz. b 50 -50x10 -9 3.0 3.0 40 -40 2.5 2.5 Drain current (µA) Photocurrent (A) IPTE [nA] ID [µA] 30 -30 2.0 2.0 1.5 1.5 20 -20 1.0 1.0 10 -10 0.5 0.5 00 00 22 44 66 88 10 10 Bias voltage (mV) SD bias [mV] Figure S7: Photocurrent and drain current measured simultaneously as a function of applied bias voltage between the source and drain contacts. We notice a linear increase of the drain current while increasing the bias, whereas the photocurrent remains constant. 20 a 140 -9 140x10 b 1000 10 -6 6 120 4 IPTE [nA] IPTE [nA] 100 100 2 Photocurrent (A) 80 -7 100 10 60 60 6 40 4 2 20 20 -8 0 0 10 10 0.0 0.5 1.0 1.5 2.0 2.5 3.0 6 8 2 4 6 8 2 4 6 0 0.5 1.0 1.5 2.0 2.5 3.0 Position (mm) 1-6 10-5 10 10 Position [mm] Pdiff [µW] Figure S8: Characterization of Device B. a) Photocurrent generated as a function of the scan position of the device. The red circles correspond to the experimental points and the red line to the Gaussian distribution fit. We obtain a beam focus FWHM of 897 µm at 2.52 THz. We observe a maximum IPTE of 136 nA for an incident power (Pin ) of 2.88 mW (Pdiff = Pin · Adiff /Afocus = 7.1 µW). Thus, the responsivity is R = IPTE /Pdiff = 19 mA/W (80.4 V/W). The sample was electrostatically doped with the left (n-doped, 60 meV) and right gate (p-doped, 50 meV) forming a pn-junction. b) Photocurrent as a function of incident power. It shows a linear trend over three orders of magnitude according to the fit displayed in dashed line. The maximum responsivity was R = 25 mA/W (105 V/W), from a photocurrent of IPTE = 145 nA, for an incident power of Pin = 2.28 mW (Pdiff = 5.6 µW, for the same beam focus at 2.5 THz and pn-junction configuration as in panel a). Using this responsivity and the Johnson √ noise that corresponds to the 4.2 kΩ measured device resistance, we find an NEP of 80 pW/ Hz. a b 40 40 -6 500x10 XPC dir. IPTE [a.u.] Y PC#1 dir. 30 30 IPTEPC[nA] 400 400 σ 3 x10 300 20 20 FWHM 200 200 10 10 100 w0 0 00 0.0 13.0 0.5 14.0 13.5 1.0 14.5 1.5 2.0 15.0 -40 -40 20 -20 00 20 20 40 40 -6 Position [mm] Position (mm) Position [μm] x10 Figure S9: Characterization of spot size for measurements on Device A. a) Photocurrent as a function of the scan position in x-direction (green) and y-direction (red) of the motorized stage onto which the device is mounted. The laser was tuned to 2.52 THz. The open dots represent the experimental data and the lines the fits according to Gaussian distributions. We obtain a FWHM of 351 (348) µm for x- (y−) scan direction. b) Calculated Gaussian distribution representing a diffraction-limit focus at 2.52 THz. The arrows indicate the standard deviation σ, FWHM and spot size w0 , corresponding to 60, 50 and 13.5% of the maximum, respectively. 21 a 0.10 b b 0.16 c c 0.32 Ⴆint,sim [A/W] Ⴆint,sim [A/W] Ⴆint,sim [A/W] 0.05 0.08 0.16 0.0 0.0 0.0 -0.05 -0.08 -0.16 -0.10 -0.16 -0.32 Figure S10: Thermoelectric simulations 50 , showing the internal responsivity (Rint,sim as defined in Eq. S11, where JE (r) represents the locally absorbed power) for a fixed device geometry and three different values of the interfacial heat conductivity Γcool . For all these three cases the device length and width dimensions are 10 and 6 µm respectively, EF = 50 meV (-50 meV) for n- (p-) doped region and with a contact resistance of 10 kΩ · µm. The Γcool is 7·104 , 2·104 and 0.7·104 W/m2 K for a), b) and c) respectively. As predicted from Eq. 2 in the main text, when decreasing Γcool , ∆T increases and therefore Rint,analyt . In fact, by calculating the Rint,analyt = IPTE /Pabs = ∆S∆T Pabs R ∆S = Aactive Γcool R as shown in the main text, we obtain excellent agreement between these calculations and the thermoelectric simulations values of Rint,sim for all three cases, namely Rint,analyt = 0.10, 0.19 and 0.33 A/W for a), b) and c) respectively. These results show that the responsivity increases (and the NEP decreases) with the square-root of the decrease in interfacial heat conductivity. 22 References (1) Tonouchi, M. Cutting-edge terahertz technology. Nature Photonics 2007, 1, 97–105. (2) Ferguson, B.; Zhang, X.-C. Materials for terahertz science and technology. Nature Materials 2002, 1, 26–33. (3) Lee, M.; Wanke, M. C. Searching for a solid-state terahertz technology. Science 2007, 316, 64–65. (4) Mittleman, D. Sensing with Terahertz Radiation. Springer-Verlag: Berlin 2003, 342, 614–617. (5) Appleby, R.; Wallace, H. B. Standoff detection of weapons and contraband in the 100 GHz to 1 THz region. IEEE Trans. Antennas Propag. 2007, 55, 2944–2956. (6) Federici, J. F.; Schulkin, B.; Huang, F.; Gary, D.; Barat, R.; Oliveira, F.; Zimdars, D. THz imaging and sensing for security applications - Explosives, weapons and drugs. Semiconductor Science and Technology 2005, 20, S266–S280. (7) Siegel, P. H. Terahertz technology in biology and medicine. IEEE Transactions on Microwave Theory and Techniques 2004, 52, 2438–2447. (8) Kleine-Ostmann, T.; Nagatsuma, T. A review on terahertz communications research. Journal of Infrared, Millimeter, and Terahertz Waves 2011, 32, 143–171. (9) Ma, J.; Karl, N. J.; Bretin, S.; Ducournau, G.; Mittleman, D. M. Frequency-division multi- plexer and demultiplexer for terahertz wireless links. Nature Communications 2017, 8, 1–7. (10) Mittendorff, M.; Shanshan Li,; Murphy, T. E. Graphene-Based Waveguide-Integrated Tera- hertz Modulator. ACS Photonics 2017, 4, 316–321. (11) Gentec-EO, Pyroelectric terahertz detectors. 2018, available at http://www.gentec- eo.com/products/thz-detectors/ . (12) Tydex, Golay cells terahertz detectors. 2018, available at http://www.tydexoptics.com/products. (13) Scontel, Superconducting Nanotechnology. 2018, available at http://www.scontel.ru/terahertz/ . (14) Sizov, F.; Rogalski, A. THz detectors. Progress in Quantum Electronics 2010, 34, 278–347. (15) Virginia Diodes. Schottky diodes terahertz detectors. 2018, available at http://www.vadiodes.com/en/products/detectors. (16) Low, T.; Avouris, P. Graphene plasmonics for terahertz to mid-infrared applications. ACS Nano 2014, 8, 1086–1101. (17) Li, Z. Q.; Henriksen, E. A.; Jiang, Z.; Hao, Z.; Martin, M. C.; Kim, P.; Stormer, H. L.; Basov, D. N. Dirac charge dynamics in graphene by infrared spectroscopy. Nature Physics 2008, 4, 532–535. (18) Koppens, F. H. L.; Mueller, T.; Avouris, P.; Ferrari, A. C.; Vitiello, M. S.; Polini, M. Pho- todetectors based on graphene, other two-dimensional materials and hybrid systems. Nature Nanotechnology 2014, 9, 780–793. 23 (19) Bonaccorso, F.; Sun, Z.; Hasan, T.; Ferrari, A. C. Graphene photonics and optoelectronics. Nature Photonics 2010, 4, 611–622. (20) Vicarelli, L.; Vitiello, M. S.; Coquillat, D.; Lombardo, A.; Ferrari, A. C.; Knap, W.; Polini, M.; Pellegrini, V.; Tredicucci, A. Graphene field-effect transistors as room-temperature terahertz detectors. Nature Materials 2012, 11, 865–871. (21) Mittendorff, M.; Winnerl, S.; Kamann, J.; Eroms, J.; Weiss, D.; Schneider, H.; Helm, M. Ultrafast graphene-based broadband THz detector. Applied Physics Letters 2013, 103, 21113. (22) Spirito, D.; Coquillat, D.; De Bonis, S. L.; Lombardo, A.; Bruna, M.; Ferrari, A. C.; Pelle- grini, V.; Tredicucci, A.; Knap, W.; Vitiello, M. S. High performance bilayer-graphene Tera- hertz detectors. Applied Physics Letters 2014, 104, 061111. (23) Zak, A.; Andersson, M. A.; Bauer, M.; Matukas, J.; Lisauskas, A.; Roskos, H. G.; Stake, J. Antenna-integrated 0.6 THz FET direct detectors based on CVD graphene. Nano Letters 2014, 14, 5834–5838. (24) Generalov, A. A.; Andersson, M. A.; Yang, X.; Vorobiev, A.; Stake, J. A 400-GHz Graphene FET Detector. IEEE Transactions on Terahertz Science and Technology 2017, 7, 614–616. (25) Bandurin, D. A.; Gayduchenko, I.; Cao, Y.; Moskotin, M.; Principi, A.; Grigorieva, I. V.; Goltsman, G.; Fedorov, G.; Svintsov, D. Dual origin of room temperature sub-terahertz pho- toresponse in graphene field effect transistors. Applied Physics Letters 2018, 112, 141101. (26) Auton, G.; But, D. B.; Zhang, J.; Hill, E.; Coquillat, D.; Consejo, C.; Nouvel, P.; Knap, W.; Varani, L.; Teppe, F.; Torres, J.; Song, A. Terahertz Detection and Imaging Using Graphene Ballistic Rectifiers. Nano Letters 2017, 17, 7015–7020. (27) Skoblin, G.; Sun, J.; Yurgens, A. Graphene bolometer with thermoelectric readout and capac- itive coupling to an antenna. Applied Physics Letters 2018, 112, 063501. (28) Mics, Z.; Tielrooij, K.-J.; Parvez, K.; Jensen, S.; Ivanov, I.; Feng, X.; Müllen, K.; Bonn, M.; Turchinovich, D. Thermodynamic picture of ultrafast charge transport in graphene. Nature Communications 2015, 6, 7655. (29) Song, J. C. W.; Rudner, M. S.; Marcus, C. M.; Levitov, L. S. Hot carrier transport and photocurrent response in graphene. Nano Letters 2011, 11, 4688–4692. (30) Lemme, M. C.; Koppens, F. H.; Falk, A. L.; Rudner, M. S.; Park, H.; Levitov, L. S.; Mar- cus, C. M. Gate-activated photoresponse in a graphene p-n junction. Nano Letters 2011, 11, 4134–4137. (31) Gabor, N. M.; Song, J. C. W.; Ma, Q.; Nair, N. L.; Taychatanapat, T.; Watanabe, K.; Taniguchi, T.; Levitov, L. S.; Jarillo-Herrero, P. Hot Carrier-Assisted Intrinsic Photoresponse in Graphene. Science 2011, 334, 648–652. (32) Cai, X.; Sushkov, A. B.; Suess, R. J.; Jadidi, M. M.; Jenkins, G. S.; Nyakiti, L. O.; Myers- Ward, R. L.; Li, S.; Yan, J.; Gaskill, D. K.; Murphy, T. E.; Drew, H. D.; Fuhrer, M. S. Sensitive room-temperature terahertz detection via the photothermoelectric effect in graphene. Nature Nanotechnology 2014, 9, 814–819. 24 (33) Guo, W.; Wang, L.; Chen, X.; Liu, C.; Tang, W.; Guo, C.; Wang, J.; Lu, W. Graphene-based broadband terahertz detector integrated with a square-spiral antenna. Optics Letters 2018, 43, 1647. (34) Goossens, S.; Navickaite, G.; Monasterio, C.; Gupta, S.; Piqueras, J. J.; Pérez, R.; Burwell, G.; Nikitskiy, I.; Lasanta, T.; Galán, T.; et al., Broadband image sensor array based on graphene- CMOS integration. Nature Photonics 2017, 11, 366–371. (35) Tielrooij, K. J.; Piatkowski, L.; Massicotte, M.; Woessner, A.; Ma, Q.; Lee, Y.; Myhro, K. S.; Lau, C. N.; Jarillo-Herrero, P.; van Hulst, N. F.; Koppens, F. H. L. Generation of photovoltage in graphene on a femtosecond timescale through efficient carrier heating. Nature Nanotechnol- ogy 2015, 10, 437–443. (36) Zuev, Y.; Chang, W.; Kim, P. Thermoelectric and Magnetothermoelectric Transport Mea- surements of Graphene. Phys. Rev. Lett. 2009, 102, 96807. (37) Tielrooij, K.-J.; Hesp, N. C. H.; Principi, A.; Lundeberg, M. B.; Pogna, E. A. A.; Banszerus, L.; Mics, Z.; Massicotte, M.; Schmidt, P.; Davydovskaya, D.; et al., Out-of-plane heat transfer in van der Waals stacks through electronhyperbolic phonon coupling. Nature Nanotechnology 2018, 13, 41–46. (38) Bistritzer, R.; MacDonald, A. H. Electronic cooling in graphene. Physical Review Letters 2009, 102, 13–16. (39) Wang, L.; Meric, I.; Huang, P. Y.; Gao, Q.; Gao, Y.; Tran, H.; Taniguchi, T.; Watanabe, K.; Campos, L. M.; Muller, D. A.; Guo, J.; Kim, P.; Hone, J.; Shepard, K. L.; Dean, C. R. One-dimensional electrical contact to a two-dimensional material. Science 2013, 342, 614–7. (40) Pizzocchero, F.; Gammelgaard, L.; Jessen, B. S.; Caridad, J. M.; Wang, L.; Hone, J.; Bøggild, P.; Booth, T. J. The hot pick-up technique for batch assembly of van der Waals heterostructures. Nature Communications 2016, 7, 11894. (41) Hwang, E. H.; Rossi, E.; Das Sarma, S. Theory of thermopower in two-dimensional graphene. Physical Review B - Condensed Matter and Materials Physics 2009, 80, 1–5. (42) Mak, K. F.; Ju, L.; Wang, F.; Heinz, T. F. Optical spectroscopy of graphene: From the far infrared to the ultraviolet. Solid State Communications 2012, 152, 1341–1349. (43) Dawlaty, J. M.; Shivaraman, S.; Strait, J.; George, P.; Chandrashekhar, M.; Rana, F.; Spencer, M. G.; Veksler, D.; Chen, Y. Measurement of the optical absorption spectra of epitaxial graphene from terahertz to visible. Applied Physics Letters 2008, 93, 131905. (44) Lundeberg, M. B.; Gao, Y.; Woessner, A.; Tan, C.; Alonso-González, P.; Watanabe, K.; Taniguchi, T.; Hone, J.; Hillenbrand, R.; Koppens, F. H. L. Thermoelectric detection and imaging of propagating graphene plasmons. Nature Materials 2017, 16, 204–207. (45) Landau, L. D.; Bell, J.; Kearsley, M.; Pitaevskiii, L.; Lifshitz, E.; Sykes, J. Electrodynamics of continuous media. Elsevier 1984, 8 . (46) Cai, W.; Moore, A. L.; Zhu, Y.; Li, X.; Chen, S.; Shi, L.; Ruoff, R. S. Thermal transport in suspended and supported monolayer graphene grown by chemical vapor deposition. Nano Letters 2010, 10, 1645–1651. 25 (47) Koh, Y. K.; Bae, M. H.; Cahill, D. G.; Pop, E. Heat conduction across monolayer and few-layer graphenes. Nano Letters 2010, 10, 4363–4368. (48) Tan, Y. W.; Kim, P. Measurement of Scattering Rate and Minimum Conductivity in Graphene. Phys. Rev. Lett. 2007, 99, 246803. (49) Lee, J. E.; Ahn, G.; Shim, J.; Lee, Y. S.; Ryu, S. Optical separation of mechanical strain from charge doping in graphene. Nature Communications 2012, 3, 1024–1028. (50) Lundeberg, M. B.; Koppens, F. H. L. Thermodynamic reciprocity in scanning photocurrent maps. 2018, In preparation. 26