Optomechanics for thermal characterization of suspended graphene

TL;DR

This paper uses optomechanical techniques to measure the thermal response of suspended graphene membranes, revealing discrepancies with classical models and proposing a boundary resistance model for better understanding.

Contribution

It introduces a noninvasive optomechanical method to characterize thermal properties of suspended graphene and accounts for boundary resistance to explain observed delays.

Findings

Measured delay times are larger than predicted by classical models.

A boundary resistance model explains the discrepancy in thermal response.

Provides a new approach for thermal characterization of atomically thin membranes.

Abstract

Thermal properties of suspended single-layer graphene membranes are investigated by characterization of their mechanical motion in response to a high-frequency modulated laser. A characteristic delay time $\tau$ between the optical intensity and mechanical motion is observed, which is attributed to the time required to raise the temperature of the membrane. We find, however, that the measured time constants are significantly larger than the predicted ones based on values of the specific heat and thermal conductivity. In order to explain the discrepancy between measured and modeled tau, a model is proposed that takes a thermal boundary resistance at the edge of the graphene drum into account. The measurements provide a noninvasive way to characterize thermal properties of suspended atomically thin membranes, providing information that can be hard to obtain by other means.

Figures

Equations

$$\frac{\mathrm{d}\Delta T}{\mathrm{d}t}+\frac{1}{\tau }\Delta T=\frac{P_{ac}}{C}{e}^{i\omega t},$$

$$z_{\omega }{e}^{i\omega t}=\alpha R{P}_{ac}\frac{e^{i\omega t}}{i\omega \tau +1}.$$

$$\tau =\frac{a^2\rho c_p}{2k},$$

$$Q_B=\frac{T_{\mathrm{sus}}-T_{\sup }}{R_B}\equiv G_B(T_{\mathrm{sus}}-T_{\sup }),$$

$$C=c_p\rho h_g\pi a^2,$$

$$R=(G_B{h}_{g}2\pi a{)}^{-1},$$

$$\tau =\frac{\rho c_{p}a}{2G_B}$$

$$\tau =\frac{\rho c_{p}a}{2G_B}=\frac{a}{2}\frac{\frac{1}{c_{1l}^2}+\frac{1}{c_{1t}^2}+\frac{\pi {\hslash }^2}{3\zeta (3)A_{uc}{k}_B^2{T}^2}}{\frac{{\bar{w}}_{1l}}{c_{1l}}+\frac{{\bar{w}}_{1t}}{c_{1t}}+\frac{\pi {\hslash }^2{\bar{w}}_{1z}{c}_{1z}}{\zeta (3)A_{uc}{k}_B^2{T}^2}}$$

$$|f_{\omega ,\mathrm{corr}}|=\frac{|f_{\omega }|}{|{\epsilon }_{\omega }|},$$

$$\angle f_{\omega ,\mathrm{corr}}=\angle f_{\omega }-\angle {\epsilon }_{\omega },$$

$$F(z(t),t)=(1+\epsilon (t))F_{\mathrm{ph}}(z(t)),$$

$$F(z(t),t)=F_{\mathrm{ph}}(z_0)+{\int }_0^t\left(\frac{\partial F_{\mathrm{ph}}}{\partial t^{\prime }}+\frac{\partial F_{\mathrm{ph}}}{\partial z}\frac{\partial z}{\partial t^{\prime }}\right)h(t-t^{\prime })\mathrm{d}{t}^{\prime }.$$

$$F(z(t),t)=F_{\mathrm{ph}}(z_0)+{\int }_0^t\frac{\partial F_{\mathrm{ph}}}{\partial z}\frac{\partial z}{\partial t^{\prime }}h(t-t^{\prime })\mathrm{d}{t}^{\prime },$$

$$F_{\mathrm{tot}}(z(t),t)=F_{\mathrm{ph},\mathrm{blue}}(z_0)+F_{\mathrm{ph},\mathrm{red}}(z_0)+$$

$${\int }_0^t\left(\frac{\partial F_{\mathrm{ph},\mathrm{blue}}}{\partial t^{\prime }}+\frac{\partial F_{\mathrm{ph},\mathrm{blue}}}{\partial z}\frac{\partial z}{\partial t^{\prime }}+\frac{\partial F_{\mathrm{ph},\mathrm{red}}}{\partial z}\frac{\partial z}{\partial t^{\prime }}\right)h(t-t^{\prime })\mathrm{d}{t}^{\prime }.$$

$$m\ddot{z}(t)+m\zeta \dot{z}(t)+Kz(t)=$$

$$F(z_0)+{\int }_0^t\left(\frac{\partial F_{\mathrm{ph},\mathrm{blue}}}{\partial t^{\prime }}+\frac{\partial F_{\mathrm{ph},\mathrm{blue}}}{\partial z}\frac{\partial z}{\partial t^{\prime }}+\frac{\partial F_{\mathrm{ph},\mathrm{red}}}{\partial z}\frac{\partial z}{\partial t^{\prime }}\right)h(t-t^{\prime })\mathrm{d}{t}^{\prime },$$

$$-{\omega }^{2}m{z}_{\omega }+i\omega m\zeta z_{\omega }+K{z}_{\omega }=i\omega F_{\mathrm{ph},\mathrm{blue}}{h}_{\omega }+\nabla F_{\mathrm{blue}}i\omega z_{\omega }{h}_{\omega }+\nabla F_{\mathrm{red}}i\omega z_{\omega }{h}_{\omega }$$

$$h(t)=1-{\mathrm{e}}^{-t/\tau },$$

$$h_{\omega }=\frac{1}{i\omega (1+i\omega \tau )}.$$

$$-{\omega }^{2}m{z}_{\omega }+i\omega m\zeta z_{\omega }+K{z}_{\omega }=\frac{F_{\mathrm{ph},\mathrm{blue}}}{1+i\omega \tau }+\frac{\nabla F_{\mathrm{blue}}}{1+i\omega \tau }{z}_{\omega }+\frac{\nabla F_{\mathrm{red}}}{1+i\omega \tau }{z}_{\omega }$$

$$-{\omega }^{2}m{z}_{\omega }+i\omega \left[m\zeta +\frac{\tau \nabla F}{1+{\omega }^2{\tau }^2}\right]z_{\omega }+\left[K-\frac{\nabla F}{1+{\omega }^2{\tau }^2}\right]z_{\omega }=\frac{A_{\mathrm{ph},\mathrm{blue}}{\epsilon }_{\omega }}{1+i\omega \tau },$$

$$z_{\omega }=\frac{A_{\mathrm{ph},\mathrm{blue}}{\epsilon }_{\omega }}{(-m{\omega }^2+i\omega m\zeta +k)(1+i\omega \tau )-\nabla F},$$

$$\mathcal{R}(z_{\omega })=A_{\mathrm{ph},\mathrm{blue}}{\epsilon }_{\omega }\frac{-m{\omega }^2-m\zeta {\omega }^2\tau +K-\nabla F}{(-m{\omega }^2-m\zeta {\omega }^2\tau +k\nabla F{)}^2+(\tau \omega (k-m{\omega }^2)+m\zeta \omega {)}^2},$$

$$\mathcal{I}(z_{\omega })=A_{\mathrm{ph},\mathrm{blue}}{\epsilon }_{\omega }\frac{-\tau \omega (k-m{\omega }^2)-m\zeta \omega }{(-m{\omega }^2-m\zeta {\omega }^2\tau +k\nabla F{)}^2+(\tau \omega (k-m{\omega }^2)+m\zeta \omega {)}^2}.$$

$$R_B=\frac{A\Delta T}{Q},$$

$$Q_{1\to 2}=U\nu \bar{w},$$

$$Q_{1\to 2}=\sum _j{U}_j{\nu }_j{\bar{w}}_j,$$

$$Q_{1\to 2}-Q_{2\to 1}=\sum _j{U}_{1j}{\nu }_{1j}{\bar{w}}_{1j}-\sum _j{U}_{2j}{\nu }_{2j}{\bar{w}}_{2j},$$

$$⟨n(\omega )⟩=\frac{1}{{\mathrm{e}}^{\hslash \omega /k_{B}T}-1},$$

$$⟨E(\omega )⟩=\hslash \omega ⟨n⟩=\frac{\hslash \omega }{{\mathrm{e}}^{\hslash \omega /k_{B}T}-1}.$$

$$\mathrm{d}{N}_{ij}=D_{ij}(\omega )\mathrm{d}\omega$$

Full text

Graphene is a 2-dimensional material with a honeycomb lattice consisting of carbon atoms Geim and Novoselov (2007). Amongst its many unusual properties, its thermal conductance has attracted major attention Pop et al. (2012); Xu et al. (2016). Extremely high thermal conductivities have been demonstrated up to 5000 W/(m$\cdot$K), well exceeding the thermal conductivity of graphite Balandin et al. (2008); Ghosh et al. (2008). These measurements were performed by Raman spectroscopy, that uses the temperature dependence of the phonon frequency Calizo et al. (2007). By measuring the thermal resistance $R$, which is the local temperature increase $\Delta T$ per unit of heat flux $\Delta Q$, one can employ analytical models of the heat transport to extract the thermal conductivity of graphene $k$. This method allowed demonstration that the thermal conductivity decreases when the number of graphene layers is increased from 2 to 4 Ghosh et al. (2010). The method has been subsequently improved, for example by better calibration of absorbed laser power Cai et al. (2010) or removing parallel conduction paths through the air Chen et al. (2010). Also the amplitude ratio between Stokes and anti-Stokes signals has been exploited Faugeras et al. (2010) as an alternative to the shift in phonon frequency. As an alternative to Raman measurements, electrical heaters Xu et al. (2014a), pump probe methods Schwarz (2016); Cabrera et al. (2015), scanning thermal microscopy Yoon et al. (2014) and temperature sensors Seol et al. (2010) have been used to study heat transport in graphene, demonstrating length dependence of the thermal conductivity Xu et al. (2014a) and a reduced thermal conductivity when graphene is supported on silicon dioxide rather than freely suspended Seol et al. (2010). Different groups have demonstrated a large variety in thermal conductivity of graphene between 600 to 5000 W/(m$\cdot$K) experimentally Balandin et al. (2008); Ghosh et al. (2008); Cai et al. (2010); Chen et al. (2010); Faugeras et al. (2010); Xu et al. (2014a); Lee et al. (2011); Dorgan et al. (2013); Chen et al. (2012); Li et al. (2014); Nika and Balandin (2016) and between 100 to 8000 W/(m$\cdot$K) theoretically Nika and Balandin (2016) making the thermal conductance of graphene a debated subject.

Besides these steady-state studies of the thermal properties of graphene, it is of interest to study its time-dependent thermal properties. This requires measurement of small temperature fluctuations in suspended graphene at frequencies in the MHz range. However, since suspended integration of temperature sensors poses problems and optical techniques for temperature measurement in suspended graphene like Raman do not offer the temperature resolution and frequency bandwidth, direct high-frequency temperature measurement in suspended graphene is difficult.

In this work, it is therefore proposed to use the thermomechanical response of suspended graphene to characterize its thermal properties at MHz frequencies. It is found that the mechanical motion is delayed by a characteristic thermal time constant $\tau$ with respect to the intensity-modulation of the laser that opto-thermally actuates the membrane. This is attributed to the time necessary for heat to diffuse through the system. The optomechanics thus provides a tool for studying the dynamic thermal properties of 2D materials. Interestingly, it is found that the measured values of $\tau$ are much higher than those expected based on literature values for the thermal conductivity $k$, specific heat $c_{p}$ and density $\rho$ of graphene. Models and measurements of drums of different diameters and on different substrates are analyzed in order to account for the large value of $\tau$. It is found that the long characteristic time is best explained by a large thermal boundary resistance at the edge of the drum.

Single-layer graphene resonators are fabricated on top of 300 nm deep dumbbell-shaped cavities (see Fig. 1). Two substrates are used, one with the cavities etched in a layer of silicon dioxide and the graphene directly transferred on top. The second substrate is coated with a layer of 5 nm chromium and 40 nm gold before graphene is transferred. This is done to help determine whether the thermal properties of the substrate influence the measured characteristic thermal time. Single layer graphene grown by chemical vapor deposition is transferred over both chips covered with a protective polymer. This polymer is dissolved and the sample is dried using critical point drying (CPD) with liquid carbon dioxide. The fluid forces in this process break one half of the dumbbell, creating a resonator on the other half with a venting channel that lets gas below the membrane escape when the vacuum chamber containing the sample is purged. The graphene is further characterized by Raman spectrocopy and atomic force microscopy (AFM) to confirm it is single layer and contamination levels are low (see Supporting Information S1).

The interferometric setup shown in Fig. 2 is used to actuate the membrane and detect the motion Bunch et al. (2007); Dolleman et al. (2015, 2016); Castellanos-Gomez et al. (2013); Zande et al. (2010). In this setup the samples are mounted in a vacuum chamber with optical access. Graphene’s motion is detected by cavity optomechanics using a red He-Ne laser, where the suspended membrane acts as moving mirror and the bottom of the cavity as a fixed back-mirror in a low-finesse Fabry-Perot cavity. The intensity of the blue laser is modulated and heats up the membrane, which will deflect due to thermal expansion. A vector network analyzer (VNA) measures the transmission between the modulation and the signal on the photodetector in a homodyne detection scheme. Frequencies between 100 kHz and 100 MHz can be measured in this setup. All measurements are performed at pressures lower than 0.02 $\mu$bar, reducing gas damping and heat transport through the gas. The red laser power was 1.2 mW incident on the sample and the blue laser power was at 0.36 mW with a large power modulation of 67% in all experiments.

Here we identify the potential source for time delay between the modulation of the blue laser and the mechanical response in the measurement setup. The block diagram in Fig. 3(a) identifies the elements and processes that play a role in actuation and detection of the membrane’s motion. The modulated intensity of the blue laser is absorbed in the graphene, generating a virtually instantaneous heating power (Fig. 3(b)) since photoexcited carriers in graphene lose their energy to phonons on timescales of a few picoseconds Wang et al. (2010). The generated heat will increase the temperature of the membrane and flow toward the substrate, resulting in a time-dependent temperature increase of the membrane, where the temperature is delayed with respect to the heating power (Fig. 3(c)). The temperature increase causes thermal expansion forces that deflect the membrane (Fig. 3(d)). At frequencies far below the resonance frequency the motion will be in-phase with the thermal expansion force, especially since the quality factor of the resonanator is typically higher than 100. The intensity modulation of the red laser due to interference effect that is used to detect the motion (Fig. 3(e)) can be regarded as instantaneous and will not cause a delay. The measurements are corrected for other delays, related to delays in the instruments (VNA, photodiode) and light path delays, using a calibration procedure discussed in the Supporting Information S3.

It is thus concluded that in the frequency range below the mechanical resonance, the delay between optical actuation and deflection in Fig. 3 is nearly completely due to the delay between heating power and temperature. A thermal system with a single time constant $\tau$, driven by an ac heating power $P_{ac}e^{i\omega t}$ can be described by the heat equation:

[TABLE]

where $\Delta T$ is the temperature difference with respect to the steady-state temperature, $C$ is the thermal capacitance and $\tau=RC$ is the thermal $RC$ product. At frequencies significantly below the mechanical resonance frequency, the thermal expansion induced amplitude $z=\alpha\Delta T$ is proportional to temperature by an effective thermal expansion coefficient $\alpha$. Solution of the heat equation gives:

[TABLE]

In section S6 of the Supporting Information a full derivation of the complex amplitude $z_{\omega}$, including the mechanical damping and inertia effects, is given based on derivations by Metzer et al’ Metzger et al. (2008); Metzger and Karrai (2004); Favero et al. (2007). This equation will be used to fit the experimental data, with the parameters $B=\alpha RP_{ac}$ and $\tau$.

An example of the measured magnitude and phase of the deflection for a resonator with a diameter of 5 $\mu$m on a cavity in silicon dioxide is shown in Fig. 4(a) and 4(b) respectively. In the 0.1 to 10 MHz range, the response is frequency dependent with a decrease in magnitude as the frequency increases. Also a phase delay is observed that increases as function of frequency. Note, that the measured phase at low frequencies is not 0, but 180 degrees. This is attributed to the small offset in the deflection that the graphene membrane has, in some membranes this was reversed in sign (indicated by 0 degrees phase at low frequencies) as shown in the Supporting Information S2. Figure 4(c) shows a measurement result which is split into a real and an imaginary part. The imaginary part of the amplitude $z_{\omega}$ can be fit by eq. 2, resulting in a value of characteristic delay time of $\tau=159$ ns, with a clearly observable maximum at radial frequency $\omega=1/\tau$. The real part of eq. 2, with the same $B$ and $\tau$, is shown in Fig. 4(c) showing a small offset with respect to the data. The offset is attributed to optical cross-talk from the modulated blue laser, which can reach the photodetector despite the optical isolation. The same effect causes the difference between the model and the magnitude and phase of the amplitude response (see Figs. 4(a) and 4(b)).

Figs. 5(a-f) show typical measurement results for drums with diameters from 2-7 $\mu$m. The optical cross-talk is significantly higher in the substrates coated with gold, which is visible as a larger offset in the real part. The resulting values of $\tau$ as function of diameter for both the silicon dioxide substrate and the gold substrate are plotted in Figs. 5(g,h), respectively. A trend is observed where $\tau$ increases as function of diameter. No significant correlation between $\tau$ and fundamental resonance frequency $\omega_{0}$ was found as shown in the Supporting Information S4, suggesting weak dependence of $\tau$ on strain.

The measured time constants in this work are significantly larger than expected based on the intrinsic properties of graphene. For example, Barton et al. Barton et al. (2012) use an expression that estimates the time constant based on the thermal properties of graphene:

[TABLE]

where $a$ is the membrane radius, $\rho$ the density of graphene, $c_{p}$ specific heat and $k$ the thermal conductivity. Using approximate values $c_{p}=600$ J/(kg$\cdot$K) (calculated in the Supporting Information S7), $k=2500$ W/(m$\cdot$K) and $\rho=2300$ kg/m3 we obtain $\tau=0.3$ ns for a 2 micron drum and $\tau=2$ ns for a 5 micron drum. The observed values of $\tau$ range between 25 to 250 ns, which is one to two orders of magnitude larger than those predicted by eq. 3. Even if the most extreme values for $c_{p}$ and $k$ are used, eq. 3 gives a lower $\tau$ than measured. The theoretical limit for $c_{p}$ is given by the Petit-Dulong law ($c_{p}=2100$ J/kg/K), and the lowest experimental literature value for $k$ is 600 W/(m$\cdot$K) Faugeras et al. (2010). It thus appears that eq. 3 cannot account for the experimental $\tau$.

Therefore we consider the possibility that the thermal conduction is limited by the substrate that supports the graphene resonator. In order to investigate this, we compare the results obtained on gold-coated and uncoated substrates. It is found that the $\tau$ on the different substrates are similar (Fig. 5(g-h)), despite the much higher thermal conductivity of the gold-coated substrate. It is thus concluded that substrate effects are not responsible for the observed value of $\tau$. This conclusion is consistent with finite element simulations (Supporting information S5) of the system.

It is well known that a thermal resistance can be present at the interface between two solids Pollack (1969); Stevens et al. (2007); Swartz and Pohl (1989); Peterson and Anderson (1973); Costescu et al. (2003); Lyeo and Cahill (2006). This effect is called interfacial thermal (or Kapitza) resistance and is caused by differences in the phonon velocities, which leads to scattering that limits the phonon transport across the interface. Several works have predicted interfacial resistances in graphene using molecular dynamics simulations Pei et al. (2012); Xu et al. (2014b). Between suspended and supported graphene a value of the boundary conductance of $2\times 10^{10}$ W/(K$\cdot$m2) was reported Xu et al. (2014b). Also, grain boundaries in graphene have been shown to cause an interfacial thermal resistance Isacsson et al. (2016). Below we argue that an interfacial thermal resistance between supported and suspended graphene could account for the unexpectedly long thermal delay times we measured.

The boundary resistance will cause the formation of a temperature discontinuity at the interface between suspended and supported graphene that can be modeled by Fourier’s law Stevens et al. (2007):

[TABLE]

where $Q_{B}$ is the boundary heat flux, $T_{\mathrm{sus}}$ the temperature in the suspended part of the graphene and $T_{\mathrm{sup}}$ temperature of the supported part. $R_{B}$ is the thermal boundary resistance and $G_{B}$ is the thermal boundary conductance. In order to estimate $G_{B}$ we use a thermal RC model, where the thermal time $\tau$ is given by de product of the heat capacity of suspended graphene $C$ and the thermal resistance $R$. It is assumed that $R$ is dominated by the interfacial thermal resistance $R_{B}$, such that $\tau$ becomes independent of $k$ of graphene:

[TABLE]

[TABLE]

where $h_{g}$ is the thickness of single layer graphene. Combining both expressions yields for the thermal time $\tau$:

[TABLE]

Using eq. 7 the thermal boundary conductance ($G_{B}=\rho c_{p}a/2\tau$) is derived from the measurements of $\tau$ as shown in Fig. 6(a), 6(b) and 6(c). This shows that the value of the thermal boundary conductance lies around 30 MW/(m${}^{2}\cdot$K). For the purpose of extracting $G_{B}$ the derived value $c_{p}=600$ J/(kg$\cdot$K) is used and a density $\rho=2300$ kg/m3. Equation 7 has been verified using finite element simulations that include a thermal boundary conductance, confirming the validity of neglecting the heat conductance $k$ (Supporting Information S5).

In order to relate the derived value of $G_{B}$ to the phonon transmission probability across the interface, the following expression is derived in Supporting Information S7:

[TABLE]

Here $c_{1j}$ is the velocity of the $j$-th phonon mode, $j=z$ for the flexural (ZA), $j=l$ for longitudinal (LA) and $j=t$ for the transverse (TA) mode. The number 1 corresponds to the suspended material. $\bar{w}_{1j}$ is the integrated transmission probability (the sum over each possible angle of incidence) of phonons over the interface, $k_{B}$ is Boltzmann constant, $T$ is temperature, $\hbar$ the reduced Planck’s constant and $A_{uc}$ is the area of the unit cell of graphene.

By using eq. (8), an average phonon transmission probability $\bar{w}$ is plotted in Figs. 6(d), 6(e) and 6(f) corresponding to the boundary conductances in Figs. 6(a), 6(b) and 6(c). The average phonon transmission probability is found to be $\bar{w}=0.3\pm 0.2$ %. Potential mechanisms that limit $\bar{w}$ include phonon interface scattering due to differences in phonon propagation velocities, boundary roughness Wen et al. (2009) and kinks Evans and Levine (2013) due to graphene edge adhesion Bunch et al. (2008). Further experimental and theoretical study of heat transport across the edge between suspended and supported graphene is needed to clarify the microscopic origins of these observations.

To summarize, a dynamic optomechanical method to measure transient heat transport in suspended graphene is demonstrated. The method does not require electrical contacts, which allows high-throughput characterization of arrays of devices. The method is used to characteristic the thermal time $\tau$ of many graphene membranes. It is found that $\tau$ is a function of diameter and its value is much larger than expected based on existing models. Measurements on gold-coated and uncoated silicon dioxide samples show similar results, showing that $\tau$ cannot be attributed to the substrate. A potential cause for the large values of $\tau$ is the presence of an interfacial thermal resistance between the suspended and supported graphene. From the measurements we determine that a thermal boundary conductance with values of $30\pm 20$ MW/(m${}^{2}\cdot$K) can account for the measurements, corresponding to a low phonon transmission probability on the order of 0.3%.

Acknowledgments

We thank Applied Nanolayers B.V. for supply and transfer or single-layer CVD graphene on our substrates. We also acknowledge useful discussions with Teun Klapwijk, Debadi Chakraborty, Daniel Ladiges and John Sader. Further, we thank the Dutch Technology Foundation (STW), which is part of the Netherlands Organisation for Scientific Research (NWO), and which is partly funded by the Ministry of Economic Affairs, for financially supporting this work. The research leading to these results also received funding from the European Union’s Horizon 2020 research and innovation programme under grant agreement No 649953 Graphene Flagship and this work was supported by the Netherlands Organisation for Scientific Research (NWO/OCW), as part of the Frontiers of Nanoscience program.

Supplemental Information

Supplemental Information: Optomechanics for thermal characterization of suspended graphene Robin J. Dolleman

Samer Houri Dejan Davidovikj Santiago J. Cartamil-Bueno Yaroslav M. Blanter Herre S. J. van der Zant Peter G. Steeneken

I Additional Experimental results and methods

S1: Graphene characterization

The graphene was examined for contamination and defects using atomic force microscopy (AFM) and Raman spectroscopy. Figure S1(a) shows a height profile obtained using tapping mode AFM on graphene supported on silicon dioxide. A small part of the graphene sheet is broken and folded over, allowing to measure the step heights between graphene on graphene and graphene on silicon dioxide. The tapping mode AFM can overestimate the step height when measuring 2D materials Nemes-Incze et al. (2008). This inaccuracy is illustrated by the disparity between the graphene-oxide and the graphene-graphene step height. However the folded-over part of the graphene allows us to say that no more than $\sim$ 1.9 nm of uniform polymer contamination could be present on the membrane. This value represents a worst-case scenario and it is very likely that the contamination is much lower. The ratio between the 2D and G peak of 2:1 in the Raman spectroscopy (see Fig. S1(b)) and the absence of a D-peak is consistent with high-quality single layer graphene with low contamination levels Ferrari et al. (2006).

S2: Example of measurement with reversed phase

In the main text it is mentioned that the response of the drum can show either a 180 degrees phase shift between laser power and mechanical motion at low frequencies, but also a 0 degrees phase is possible. This is equivalent to a sign change of the real and imaginary part of the amplitude. Although 180 degrees is usually observed in our measurements, occasionally a 0 phase is observed as shown in Fig. S2. Here, the sign of real and imaginary part have changed with respect to the examples in the main text. The cause of this difference in phase are small asymmetries in the initial state which will cause either upward or downward motion due to the thermal expansion.

S3: Calibration procedure

In order to correct the intrinsic phase shifts in our measurement setup, we directly point the blue laser to the photodiode to obtain an calibration curve for our system (Figure S3). This can be corrected by deconvolution of the measured response with this calibration curve, which is done by expressing the blue laser modulation parameters as a frequency-dependent phasor $\varepsilon_{\omega}$. Since the calibration was taken at discrete frequencies, a cubic interpolation was used to make sure the frequencies match the ones from the measurement that needs to be corrected. Now one can deconvolve the measured frequency response function $f_{\omega}$ using:

[TABLE]

[TABLE]

where $f_{\omega,\mathrm{corr}}$ is the corrected frequency response function of our measurement.

S4: Correlation between resonance frequency and $\tau$

Figure S4 shows scatter plots between the resonance frequency and characteristic thermal time $\tau$ extracted from the measurements. Since both the resonance frequency and characteristic thermal time are correlated to diameter, correlations between the two variables should only be determined for the same diameter. We found low correlations close to zero with outliers at -0.24 for the 6 micron diameter drums on gold and 0.14 for the 5 micron diameter drums on gold. The low correlations and the low agreement between different diameters suggests that transient thermal transport is not strongly related to the strain present in the graphene resonators.

II Additional modeling results

S5: Finite element simulations of graphene on a silicon dioxide substrate

In order to examine the impact of the silicon dioxide substrate on the heat transport, we use COMSOL Multiphysics to model the graphene on top of the cavity and estimate the delay function $h(t)$. A simulation result of $h(t)/P_{\mathrm{laser}}$, where $P_{\mathrm{laser}}$ is the incident laser power (assuming $2.3\%$ absorption of optical power), can be seen in Fig. S5. This simulation predicts that the heat transport is more complex than expected. A very fast increase in temperature is observed with a time constant that is in the order of 0.5 ns. This is followed by a much slower exponential increase in temperature, which can be fitted with a single exponential to obtain a time constant of 32.6 ns. The fast time constant should not be observed in our measurement, since the cut-off frequency $2\pi\omega_{c,\mathrm{fast}}=1/\tau_{\mathrm{fast}}\approx 320$ MHz is much larger than the bandwidth in our measurements. The slow time constant can be observed in our measurement, since the cut-off is in a measurable frequency range and lower than the resonance frequency. This could be the thermal relaxation time found in the measurements.

Figure S6 shows $h(t)$ for different material parameters. From fig. S6(a) we conclude that changes in the thermal conductivity or specific heat of the graphene membrane affect the fast time constant $\tau_{\mathrm{fast}}$, the slow time constant $\tau_{\mathrm{slow}}$ remains unchanged. We conclude that the observed time constant $\tau$ in our measurements does not depend on the properties of graphene itself. For the thermal contact resistance between the graphene and silicon dioxide interface shown in fig. S6(b), we draw the same conclusion.

Figure S6(c) shows a different situation; if the thermal properties of the silicon dioxide layer are changed, the fast time constant remains unchanged. However, the slow time constant changes a lot, from 32 ns to 66 ns if the thermal conductivity is changed from 1.4 W/m/K to 0.7 W/m/K. We therefore conclude that the slow time constant in the simulation depends only on the properties of the substrate. Due to the poor thermal properties of this layer, it takes much longer to reach thermal equilibrium than expected if only the graphene itself is considered. Figure S6(d) shows the diameter dependence on $h(t)$, again the slow time constant is hardly affected, therefore this model does not account for the diameter dependence in our measurement. Since the slow time constant modeled here only depends on the properties of the substrate, it should significantly change value when the silicon dioxide is replaced with gold. Since the values simulated here are not diameter dependent, this model shows that the experimental observations cannot be explained by the thermal properties of the substrate.

Figure S7 shows a simulation with identical parameters as in Fig. S5 with the addition of a limited thermal boundary conductance $G_{B}=$ 40 MW/(m${}^{2}\cdot$ K). The interfacial thermal resistance dominates the heat transport in this situation as illustrated by the uniform temperature in the drum. The materials outside the drum do not raise in temperature significantly. This validates the simple model in the main part of this work, where only the boundary conductance is considered.

S6: Derivation of expression for motion of drum actuated by modulated laser

The derivations below follow the methods used by Metzger et al.Metzger et al. (2008); Metzger and Karrai (2004); Favero et al. (2007) closely. The derivation is repeated here to show that the use of a separate red laser to read out the motion does not affect the amplitude response. We assume there is a weakly intensity-modulated laser actuating the graphene membrane, the force induced by this laser $F$ can be written as:

[TABLE]

where $F_{\mathrm{ph}}$ is a photo-induced force (which can be photo-thermal, radiation pressure or radiometric pressure) that is exerted on the compliant graphene membrane and $\varepsilon$ a modulation parameter assumed to be much smaller than 1. In general, the graphene membrane will show a delayed response in its deflection due to this force, but responds with a certain time delay $\tau$. This time-delay can be described by a function $h(t-t^{\prime})$ that leads to the following formulation for the force:

[TABLE]

In the case of constant illumination, this equation becomes:

[TABLE]

which is valid for the red laser used in the interferometric detection. The combined action of the red and blue laser gives for the total force:

[TABLE]

The equation of motion that needs to be solved now reads:

[TABLE]

where $z(t)$ is the amplitude on a generalized coordinate on the membrane, which we assume is the amplitude detected by the red laser; $m$ is the effective mass, $\zeta$ the effective damping coefficient and $K$ the effective stiffness. We drop the term $F(z_{0})$, since this leads only to a static deflection with no time-dependence. We assume that the amplitude is very small and that the terms $\partial F_{\mathrm{ph}}/\partial z$ can be approximated by a constant $\nabla F$. Using the properties of Laplace transforms for convolutions we can now write eq. II in the frequency domain:

[TABLE]

We assume the shape of the delay function is of the exponential type:

[TABLE]

which has the Laplace transform:

[TABLE]

Inserting this into the equation of motion gives:

[TABLE]

. Now we split the effective actuation force into $F_{\mathrm{ph,blue}}=A_{\mathrm{ph,blue}}\varepsilon_{\omega}$ where $A_{\mathrm{ph,blue}}$ represents the force from the DC blue laser power and $\varepsilon_{\omega}$ is the modulation parameter. $\varepsilon_{\omega}$ is made constant as function of frequency by the calibration method described in section S3. Regrouping in terms of omega gives:

[TABLE]

where we combined $\nabla F_{\mathrm{blue}}$ and $\nabla F_{\mathrm{red}}$ into $\nabla F=\nabla F_{\mathrm{blue}}+\nabla F_{\mathrm{red}}$. The solution for the amplitude is:

[TABLE]

with real and imaginary part:

[TABLE]

In the limit where $1/\tau\ll\sqrt{k/m}$ the imaginary part shows a local maximum at a radial frequency of approximately $1/\tau$. In this case, the frequency response function is $\omega\tau/(1+\omega^{2}\tau^{2})$, which is used in this work to fit to the imaginary part of the measured response.

Derivations of interfacial thermal resistance, specific heat and characteristic thermal time for 2D materials

S7: Derivation of expression for $\tau$ (eq. 8 in main text)

The interfacial thermal resistance $R_{B}$ can be determined by using the heat flux:

[TABLE]

where $A$ is cross sectional area of the boundary, $\Delta T$ the temperature difference and $Q$ the heat flux. The first step in determining the interfacial resistance, is thus to determine the heat flux that crosses the interface. The heat flux that crosses from interface 1 to interface 2 can be expressed by Gray (1981); Peterson and Anderson (1973):

[TABLE]

where $U$ is the total energy per unit volume of the heat carriers, $\nu$ the velocity at which they propagate and $\bar{w}$ the probability that the heat carriers transmit over the interface. In the calculation of thermal interfacial resistance, the difficulty lies in calculating the transmission probability $\bar{w}$, while the calculation of energy and propagation velocity is quite straigtforward.

Our approach is thus, to calculate the energy and velocity and use that to estimate the value of $\bar{w}$ from the measurement. In order to do this, it is assumed that the heat in graphene is carried by phonons and that all the heat is carrier by three acoustic phonon polarizations, the longitudinal (LA), transverse (TA) and flexural (ZA). The LA and TA branch are far below the Debye temperature of 2100 K due to their large velocities, but the ZA branch will be fully thermalized since its Debye temperature is at 50 K Hone (2001). The contribution to the heat flux of each polarization can be added to obtain:

[TABLE]

and the total heat flux becomes:

[TABLE]

the index $ij$ now describes the material ($i=1$ for suspended, $i=2$ for supported graphene) and phonon mode $j$.

To calculate the energy $U_{ij}$, we start from the Bose-Einstein distribution to find the average phonon number $\langle n(\omega)\rangle$ at a fixed frequency $\omega$:

[TABLE]

where $\hbar$ is the reduced Plack’s constant, $k_{B}$ the Boltzmann constant and $T$ temperature. The energy carried by each phonon is $\hbar\omega$, therefore we can write the average energy $\langle E(\omega)\rangle$ that phonons have at this frequency:

[TABLE]

The number of states (density of states $D(\omega)$) that is accessible to the system between the frequencies $\omega$ and $\omega+\mathrm{d}\omega$ is defined by:

[TABLE]

and this makes the total energy $V$ in the system:

[TABLE]

Next, one has to know for the density of states, how many modes are available in the momentum space. In a 2-dimensional crystal, if we know the size of the system $A$, the uncertainty in momentum is $(2\pi\hbar)^{2}/A$ and the number of modes available to the system becomes:

[TABLE]

where one can integrate over circles with circumference $2\pi k$ to obtain:

[TABLE]

Using $\mathrm{d}N=D(\omega)\mathrm{d}\omega$ we obtain for the total energy in the system:

[TABLE]

which is divided by the total volume of the system to obtain for $U_{ij}$:

[TABLE]

to perform the integration, it is necessary to use the dispersion relation that relates the frequency to the wavenumber. Since the flexural phonons have different properties than the transverse and longitudinal phonon, these will have to be analyzed separately in the sections below.

II.0.1 Longitudinal and tranverse modes

For the longitudinal and transverse acoustic phonons, we can write the linear dispersion relationship:

[TABLE]

where $c_{ij}$ is the propagation velocity of the phonons, since $\nu=\frac{\mathrm{d}\omega}{\mathrm{d}k}=c_{ij}$. Substitution into eq. S27 gives:

[TABLE]

Using this, we can write for the heat flux over an interface of area $A$ from material 1 to material 2:

[TABLE]

To solve the frequency integral, we assume $k_{B}T\ll\hbar\omega_{D}$ and using a coordinate transform $x=\hbar\omega/k_{B}T$:

[TABLE]

We assume that the transmission probability frequency-independent. This results in:

[TABLE]

and for the heat flux from material 2 to material 1:

[TABLE]

II.0.2 Flexural mode

When strain is present in a 2D lattice, the dispersion of the flexural phonons can be written as Lifshitz (1952):

[TABLE]

which has 4 solutions for $k$, however implying the conditions $s_{iz}>0$, $c_{iz}>0$, $k>0$ and enforcing that $k$ must be a real number, we only have one solution:

[TABLE]

this can be substituted in eq. S27 to obtain:

[TABLE]

Note, that this equation converges either to the expression for linear dispersion if $c_{iz}^{4}\gg 4s_{iz}\omega^{2}$ or to the expression for quadratic dispersion if $c_{iz}^{4}\ll 4s_{iz}^{2}\omega^{2}$. The propagation velocity becomes:

[TABLE]

substituting eq. S35 gives:

[TABLE]

which in the limit case of purely quadratic dispersion $c_{iz}=0$ becomes $\nu_{iz}=2\sqrt{\omega s_{iz}}$ and in the case of linear dispersion $s_{iz}=0$ becomes $\nu_{iz}=c_{iz}$.

The analysis can be simplified by assuming either high or low strains, which should follow from our experiments. It can be seen, that the condition:

[TABLE]

allows us to use descibe the heat transport of the ZA branch by a quadratic dispersion without strain, while the condition:

[TABLE]

allows one to use a linear dispersion for the ZA branch. From the resonance frequencies in our experiments, we can estimate the strain present in the drum resonators:

[TABLE]

assuming $\rho h=$7.7\text{\times}{10}^{-7} kg/m2 and $Eh=$ 340 N/m, we find the lowest observed strain $\epsilon_{\mathrm{low}}=$1.026\text{\times}{10}^{-5} and the highest observed strain $\epsilon_{\mathrm{high}}=$1.71\text{\times}{10}^{-4}$$. Expressions for coefficients $s_{iz}$ and $c_{iz}$ are given by Lifshitz Lifshitz (1952):

[TABLE]

where $\kappa$ is the bending rigidity, which is $\kappa=$1\text{\times}{10}^{-19}$$ J for single-layer graphene, and:

[TABLE]

where $u$ is the dilatation and $\lambda$, $\mu$ are the Lame parameters. Now we can calculate the coefficients for the lowest strain:

[TABLE]

[TABLE]

from which we conclude that for each drum measured in this work the condition in eq. S40 holds. Due to the low velocities the Debye temperature of the flexural phonons is much lower than the in-plane phonons. Therefore, we have to write the heat flux as:

[TABLE]

where $\theta$ is the Debye temperature, since we are above the Debye temperature $e^{x}\approx 1+x$:

[TABLE]

the total number of states in the system is:

[TABLE]

the Debye temperature becomes:

[TABLE]

here $N/\pi a^{2}$ is the number of states per unit square, which is limited by the area of the unit cell $A_{uc}=$5\text{\times}{10}^{-20}$$ m2. The heat flux from the ZA mode now becomes:

[TABLE]

and the total heat flux now becomes:

[TABLE]

The condition $Q_{tot}=0$ has to apply if $T_{1}=T_{2}$, this implies that:

[TABLE]

which makes the heat flux:

[TABLE]

This can be linearized for small temperature differences $\Delta T$ to obtain:

[TABLE]

and the boundary resistance is directly obtained from eq. S16:

[TABLE]

II.0.3 Specific heat

We can also calculate the specific heat $c_{p}$ by starting from equation S29:

[TABLE]

[TABLE]

which is valid for the LA and TA branches. For the ZA phonons we have to take the high temperature limit:

[TABLE]

using $\mathrm{e}^{x}\approx 1+x$:

[TABLE]

now by substituting:

[TABLE]

the energy density becomes:

[TABLE]

[TABLE]

Now we find:

[TABLE]

[TABLE]

In the main text we found the expression:

[TABLE]

Substituting eqs. S55 and S64 gives:

[TABLE]

this result was used in the main text to estimate the average phonon transmission probability $\bar{w}$. This average is defined as the situation where all branches have equal transmission probability: $\bar{w}=\bar{w}_{1t}=\bar{w}_{1l}=\bar{w}_{1z}$.