A graphite thermal Tesla valve driven by hydrodynamic phonon transport | Nature

Nature (2024)Cite this article

Metrics details

The Tesla valve benefits the rectification of fluid flow in microfluidic systems1,2,3,4,5,6 and inspires researchers to design modern solid-state electronic and thermal rectifiers referring to fluid-rectification mechanisms in a liquid-state context. In contrast to the rectification of fluids in microfluidic channels, the rectification of thermal phonons in micro-solid channels presents increased complexity owing to the lack of momentum-conserving collisions between phonons and the infrequent occurrence of liquid-like phonon flows. Recently, investigations and revelations of phonon hydrodynamics in graphitic materials7,8,9,10 have opened up new avenues for achieving thermal rectification. Here we demonstrate a phonon hydrodynamics approach to realize the rectification of heat conduction in isotopically enriched graphite crystals. We design a micrometre-scale Tesla valve within 90-nm-thick graphite and experimentally observe a discernible 15.2% difference in thermal conductivity between opposite directions at 45 K. This work marks an important step towards using collective phonon behaviour for thermal management in microscale and nanoscale electronic devices, paving the way for thermal rectification in solids.

In microfluidics, the Tesla valve, invented and patented by Nikola Tesla in 1920 (ref. 1), is widely used to manipulate fluid flow2,3,4,5,6. The inherent asymmetry of the Tesla valve leads to the distinct fluid-flow behaviours in opposite directions. In the reverse flow direction, the viscous fluid encounters notable impedance, whereas in the forward direction, the flow experiences considerably less resistance. Thus, it realizes the rectification of fluid flows without the need for any moving parts2,6.

Analogous to macroscopic phenomena in fluid dynamics, the collective motion of electrons or phonons manifests extraordinary hydrodynamic electrical11,12,13 or thermal14,15,16 conduction in solids. In particular, phonon hydrodynamics shows peculiar characteristics of collective behaviours, such as the second sound8,9,17, phonon Poiseuille flow10,18,19, phonon vortex20,21,22 and phonon Knudsen minimum23,24,25,26. Despite the recent observations of hydrodynamic thermal phenomena, contemporary phonon hydrodynamics remains in an exploratory phase, with practical applications yet to be realized. This is primarily attributed to stringent observation conditions7,27 and the limited availability of materials capable of hydrodynamic phonon transport28,29.

On the other hand, the rectification of thermal conduction has promising application in engineering thermal energy, demonstrating potential capability in the design of thermal diodes30,31,32. The thermal-rectification effect has been claimed in many solid systems comprising two components with different temperature-dependent thermal properties33,34,35. It allows stronger heat flux in one direction and weaker heat flux in the opposite direction. However, induction of external treatments is required to trigger a phase transition of the material33 or atomic-scale synthesis is required to avoid further effects arising from heat conduction across interfaces between heterogeneous compounds34. Furthermore, trapezoidal and asymmetric nanoporous structures have been devised to achieve rectification of thermal conduction in homogeneous material systems36,37,38. However, either a relatively large temperature difference is required to achieve the thermal rectification36,37 or the directivity of the thermal behaviour is limited to the demonstration of local rectification effects within nanometre/submicrometre scales36,38.

In this study, we designed asymmetric structures inspired by the Tesla valve used in hydrodynamics (Fig. 1a). We fabricated micrometre-scale samples using graphite, a material known for its inherently robust hydrodynamic phonon behaviour8,10,24. Using the thermoreflectance technique, we examined the thermal properties and rectification effects exhibited by the suspended graphite Tesla valve across the temperature range of 10–300 K. Analogous to hydrodynamics, the thermal rectification observed in solid graphite is attributed to hydrodynamic-like phonon flows in the Tesla valve. The concept of a hydrodynamic phonon Tesla valve in this work holds promise for achieving thermal rectification within complex micrometre-scale structures, thereby introducing new prospects for thermal manipulation using phonon hydrodynamics.

a, Optical microscope image of graphite Tesla valves in forward (top) and reverse (bottom) directions. Scale bar, 10 μm. b, Schematic of phonon hydrodynamic flows in a graphite Tesla valve in the reverse case. Arrows indicate the paths of the phonon hydrodynamic flows. Inset, schematic of the pump–probe thermoreflectance measurement. c, Example fitting of the experimental curve of the TDTR signal by FEM simulations in a forward structure at 45 K. The thermal conductivity of the graphite Tesla valve was extracted by obtaining the best fit. d, Histogram of recorded thermal decay times during one measurement. e, Temperature distribution in the Tesla valve structure in the forward direction simulated by FEM. a.u., arbitrary units.

Source Data

In graphite, the out-of-plane phonon modes induced by the interlayer interactions of graphene are largely populated at the Γ point39,40. This leads to intense anharmonic collisions among phonons with small wavevectors, consequently promoting frequent normal scattering in a momentum-conserving manner28,41. The strength of such conservation of phonon momentum plays a crucial role in the maintenance of hydrodynamic phonon flows and corresponds closely to the viscosity in fluid flows7. On the other hand, carbon-based materials have a relatively high Debye temperature (approximately 2,000 K)42 and large-wavevector phonon modes are not populated at temperatures much below their Debye temperature. Consequently, graphite exhibits a temperature window in which the hydrodynamic phonon flow occurs owing to predominance of normal scattering8,17,24,43.

Recent works have reported that the hydrodynamic phonon transport is strongly hindered by external defect scattering, including isotopes7,10 and vacancies24. To minimize these effects, we used isotopically enriched graphite crystals (comprising 99.98% 12C and 0.02% 13C). Concurrently, electron backscatter diffraction (EBSD) analysis confirmed the single-crystalline structure of the graphite samples (Methods and Extended Data Fig. 1). The observed single crystallinity of the samples creates a more conducive environment for hosting collective phonon behaviour, minimizing further external scattering arising from polycrystalline grain boundaries.

We performed experiments based on a well-developed fabrication and measurement system tailored for micrometre structures of exfoliated materials. The main heat-flow channel of the graphite Tesla valve had length 26 μm and width 4.5 μm, with the bent channel having a width equivalent to that of the main channel. Both the forward and reverse structures were fabricated on a single graphite flake with a uniform thickness of 90 nm (Extended Data Fig. 2). A reversal in heat-flow direction was realized by inverting the heat source and the heat sink after carefully excluding the impact of different resistances at the boundaries of these two ends (details in Methods and Extended Data Fig. 3). We suspended the structures to avoid heat loss through the substrate and placed the samples in a cryostat with high vacuum (<10−4 Pa) to facilitate low-temperature measurement down to 10 K. The thermal properties of the graphite Tesla valve were measured using the microsecond time-domain thermoreflectance (μ-TDTR) method (Fig. 1b) and the values of thermal conductivity were obtained by processing the data with finite element method (FEM) simulations (Fig. 1c–e and details in Methods). We prepared three graphite samples, namely, S1–S3, and the following experimental results are taken from the last sample (S3) unless indicated otherwise.

We measured the thermal conductivity (κ) of the graphite Tesla valve for both forward and reverse directions in the temperature range 10–300 K. At 10 K, we obtained comparable values of κ for both forward (κf = 4.18 ± 0.11 Wm−1 K−1) and reverse (κr = 4.14 ± 0.04 Wm−1 K−1) directions. The negligible difference indicates that, at low temperatures, at which internal phonon–phonon scattering is absent, ballistic phonons excited by the heat source in random directions can traverse the entire structure—including both the main and the bent channels—unimpeded by anything except the boundaries of these channels.

Figure 2 illustrates the Monte Carlo simulations that were conducted to investigate the behaviour of the graphite Tesla valve in the absence of the hydrodynamic effect. Specifically, phonons were initiated from the heat source, exhibiting a randomized distribution of frequencies in accordance with Planck’s distribution at a temperature of 10 K. Their trajectories covered the entire Tesla valve channels, terminating at the cold ends, as depicted in Fig. 2a,c for the forward and reverse cases, respectively. The corresponding heat-flux maps, presented in Fig. 2b,d, manifested a flux mainly through the main channel from the hot to cold ends in both scenarios. Details of the Monte Carlo method can be found in Methods and our previous works44,45,46.

a–d, Phonon paths and heat-flux maps of the Tesla valve in forward (a,b) and reverse (c,d) directions at 10 K. Scale bar, 5 μm. e, Phonon mean-free-path spectra of the Tesla valve in both directions at 10 K. The peak positions correspond to the width and length of the channel. Inset, diodicity as a function of temperature from Monte Carlo simulations showing no rectification in the ballistic regime.

Source Data

The analysis of the phonon mean free paths spectrum within the Tesla valve structures is shown in Fig. 2e. The spectra suggest that the predominant confinement of phonon mean free paths occurs within the characteristic size, primarily corresponding to the width of the channels. Such confinement is attributed to the dominant boundary scattering observed in the ballistic regime for both forward and reverse directions. Consequently, this confinement explains the negligible difference between κf and κr at 10 K under particle phonon transport approximation. These simulations help to ensure that the given geometries create no thermal rectification outside the hydrodynamic regime.

With a gradual increase in temperature, anharmonic lattice vibrations intensify, leading to stronger phonon–phonon interactions47,48. In contrast to Umklapp scattering, which violates the conservation of phonon momentum in k-space, normal scattering preserves phonon momentum during collisions, closely resembling the collision dynamics observed between fluid molecules. Analogous to the viscous fluid flow in a finite-width channel with a parabolic flow profile, known as Poiseuille flow, the hydrodynamic flow of phonons is shaped by the interplay between normal scattering and boundary scattering at the channel walls. To establish the hydrodynamic phonon flow, or phonon Poiseuille flow, the following conditions must be met49: predominance of normal scattering to realize the Poiseuille flow of phonons within the channel under a macroscopic description (that is, λN ≪ W, in which λN is the mean free path of normal scattering and W is the channel width) and sufficiently weak resistive scattering, including Umklapp and defect scattering, to minimize phonon momentum loss (that is, λNλR ≫ W2, in which λR is the mean free path of resistive scattering).

Previous theoretical studies predicted the hydrodynamic phonon flow in graphite within both temperature (20–90 K) and width (2–20 μm) windows24,43. A recent experimental study reported the observation of Poiseuille flow of phonons in isotopically enriched graphite with channel widths of 3.3 and 5.5 μm, respectively10. Our Tesla valve channels, with a width of 4.5 μm, are designed within the hydrodynamic width window. This design ensures the hosting of hydrodynamic phonon flows within the specified temperature regime, thereby enabling the phonon hydrodynamic effects.

As the temperature increases above 20 K, in the forward case (Fig. 3a), the collective flow of phonons maintained by the normal scattering passes predominantly through the main channel of the Tesla valve. Only a smaller fraction of phonons deviates into the bent channel, entering at an angle opposite to the direction of heat flow, as evidenced by our thermal energy map simulations (Extended Data Fig. 4). However, in the case of the reverse direction (Fig. 3b), two hydrodynamic phonon flows diverge from the heat source. One passes through the main channel and the other is diverted to the bent channel with a straight entrance. The convergence of these two hydrodynamic phonon flows subsequently causes the redirection of the phonon momentum and the loss of the heat flux at the exit of the bent channel, as depicted in Fig. 1b. Such a nonlinear thermal behaviour would be established only in the hydrodynamic regime and consequently introduce a resistance to the thermal transport, resulting in a reduced thermal conductivity compared with that observed in the forward case. In Fig. 3c, an appreciable difference between κf and κr is evidenced from 25 to 60 K, the temperature range situated within the hydrodynamic temperature window24,43.

a,b, Optical microscope images of the graphite Tesla valve in forward and reverse directions, respectively. Arrows indicate the directions of heat flows. Scale bars, 5 μm. c, Thermal conductivity of the graphite Tesla valve in forward (cyan triangles) and reverse (red triangles) directions. Error bars denote the standard deviations (σ) of κ values from 3 to 10 measurements for each structure (forward or reverse) at a given temperature.

Source Data

In fluid dynamics, the efficiency of rectification (or diodicity) using a Tesla valve is expressed as the ratio of the flow resistance in the reverse direction to that in the forward direction2. Similarly, we defined the diodicity of the heat flow in our graphite Tesla valve structures as the ratio of the thermal conductivity in the forward direction to that in the reverse direction, D = κf/κr. Therefore, the rectifying effect in the graphite Tesla valve is manifested when thermal conduction in the forward direction is less resistive than that in the reverse direction, that is, when D is greater than 1.

Figure 4a shows the thermal diodicity of the graphite Tesla valve from 10 to 300 K. Below 20 K, the diodicity remains at the baseline value of 1, suggesting that phonon transport through the Tesla valve is the same in both directions. As the temperature increases, phonons exhibit a collective motion maintained by momentum-conserving normal scattering. At around 45 K, phonon hydrodynamic fluxes formed in the graphite Tesla valve channels, in which resistive Umklapp scattering is still scarce. As discussed earlier, the asymmetric structure of the Tesla valve facilitates unhindered flow of phonons in the forward direction, whereas the flow is more resistive in the reverse direction. Consequently, the thermal conductivity for the forward case is 15.2% higher than that for the reverse case, as indicated by the peak diodicity at 45 K.

a,b, Thermal diodicity (the ratio of thermal conductivity in the forward, κf, and reverse, κr, directions) as a function of temperature in graphite Tesla valve samples S1–S3 (a) and the silicon Tesla valve (b). The orange region shows the calculated phonon hydrodynamic temperature window of a 4.5-μm-wide graphite channel in our previous work43. The top schematics illustrate different phonon-transport behaviours in graphite within the hydrodynamic regime and in silicon within the diffusive regime. The error bars in the plots are calculated as \(\left({\kappa }_{{\rm{f}}}/{\kappa }_{{\rm{r}}}\right)\sqrt{{\left({\sigma }_{{\rm{f}}}/{\kappa }_{{\rm{f}}}\right)}^{2}+{\left({\sigma }_{{\rm{r}}}/{\kappa }_{{\rm{r}}}\right)}^{2}}\), in which σf and σr are the standard deviations of κf and κr, respectively.

Source Data

Further experiments were performed previously on graphite Tesla valve samples S1 and S2 with thicknesses of approximately 90 nm and 140 nm, respectively. At 40 K, we obtained the diodicity values of approximately 1.14 and 1.11 in S1 and S2, respectively (the thermal conductivities are provided in Extended Data Fig. 5), thereby corroborating the reproducibility of the thermal-rectification phenomenon within our graphite Tesla valve.

Notable thermal rectification in the 25–60-K range was independently observed in all three samples (S1–S3). This observation is consistent with the predicted hydrodynamic window of 20–65 K (ref. 43). As the temperature continues to increase, phonons populate higher frequency (or larger wavevector) states, resulting in the enhanced Umklapp scattering that breaks phonon momentum28. Thus, the weakening of the hydrodynamic phonon behaviour occurs in our experiments, as shown by the gradual disappearance of thermal rectification at temperatures above 60 K.

In the realm of classical diffusive phonon transport, heat flow obeys Fourier’s law, resulting in the negligible thermal rectification as evidenced by a constant value of diodicity around 1 up to 300 K. High-temperature experiments conducted on samples S1–S3 reproduced a similar trend, albeit with unavoidable error bars owing to the mechanical instability; however, the diodicity averages about 1 for these three graphite Tesla valve samples, indicating the absence of thermal-rectification effects in the diffusive regimes (see detailed experimental error analysis in Methods).

It is noteworthy that, as well as the asymmetric structure (Tesla valve), the asymmetry (or nonlinearity) of thermal transport in the hydrodynamic temperature regime is a crucial factor in ensuring thermal rectification. To demonstrate this fact, we extended our investigations to silicon, a material conventionally characterized as non-phonon-hydrodynamic with minimal momentum-conserving normal scattering, in which Umklapp scattering predominantly governs most of the phonon–phonon interactions in a momentum-destroying manner14,50.

In a 145-nm-thick suspended silicon thin membrane, we fabricated and measured Tesla valve structures identical to those constructed with graphite (Extended Data Fig. 6). As shown in Fig. 4b, we examined the thermal-rectification phenomenon across the entire temperature range. In contrast to the hydrodynamic-like heat flows in graphite in the hydrodynamic regime, the heat flux in silicon undergoes substantial damping as a result of the momentum-destroying effects of Umklapp scattering. Consequently, heat transports diffusively in both the forward and reverse directions and the diodicity values stay constant around 1. This was confirmed in another supplementary silicon sample at hydrodynamic temperatures (Extended Data Fig. 7). At low temperatures, at which phonons with long mean free path may cause quasi-ballistic thermal transport, we observed diodicity around 1, consistent with the previous conclusion that the rectification is hard to see without nonlinear effects51. Thus, hydrodynamic phonon behaviour is required to observe thermal rectification in a Tesla valve.

As discussed previously, the micrometre-scale design of the channels plays a pivotal role in shaping collective phonon motions and Poiseuille heat flows. For the reverse direction of our phonon Tesla valve, our design considers straight inlets for the bent channels aligned with the direction of the heat flux induced by the source. Also, 21-μm-long straight paths are designed to establish pronounced phonon hydrodynamic flows in the bent channels. A recent study explored the rectification of electron flows in graphene Tesla valve structures under the electron hydrodynamic framework52. The observed weaker rectification effect of about 5% in their study could be potentially explained by overlooking the two critical factors mentioned above, that is, the straight entrance and the length of the bent channels. Both factors are well-established requirements in fluid dynamics Tesla valves3,6 and are probably required for observing thermal rectification in solid Tesla valves using phonon hydrodynamics.

Further enhancement of thermal rectification can be achieved by connecting several segments of Tesla valve units in series, a well-established concept in fluid dynamics4. However, the current limitations in fabrication confine the actual suspended graphite structures to the flake scale of approximately 100 × 100 μm2. Future experimental attempts to improve the graphite phonon Tesla valve involve exploring wider channel widths under hydrodynamic window conditions and the extension of the exit lengths of the bent channels. These strategies are being pursued with the goal of increasing the diodicity value for improved thermal-rectification effects in the graphite Tesla valve and will be explored in our future research.

We first prepared commercial 2.5-μm-thick SiO2 on Si wafers and baked them with air plasma for the surface treatment. Graphite flakes were transferred to the SiO2/Si substrate by mechanical exfoliation using Scotch Tape. We used electron-beam lithography to pattern the Tesla valve profiles on the initial flake and deposited an Al thin film to mask them by electron beam physical vapour deposition. We then used O2 plasma to etch the exposed flake in a reactive ion etching system and the designed structures were obtained after the removal of the top Al thin film. In the centres of the two symmetrical Tesla valves, we deposited 70-nm Au pads as transducers for performing the thermoreflectance measurements in both the forward and reverse structures. Finally, the graphite Tesla valves were suspended by the Au pads deposited on both sides of the structures in combination with the underlying SiO2 etching using vapour-phase hydrofluoric acid.

We performed EBSD analysis on the entire 50 × 150-μm2 flake (Extended Data Fig. 1a) to investigate the crystallinity of our graphite samples. As shown in Extended Data Fig. 1b, the flake was irradiated with the electron beam at a tilted angle of 70° to the wafer. The electron diffraction pattern recorded the information of each crystal plane. The inverse pole maps (Extended Data Fig. 1c–e) show the consistency of crystal orientation (identical colour) for all three directions, namely, the normal direction, the transverse direction and the reference direction. As a result, we confirmed that the graphite samples used in this work are single crystalline.

In-plane quasi-ballistic thermal conduction in materials with a hexagonal lattice is reported to be dependent on the different lattice orientations of armchair and zigzag53. To exclude this effect, we fabricated our samples on the single-crystalline graphite flake so that the main and the bent channels in the forward structure are parallel (or have 180° symmetry) to their counterparts in the reverse case and identical lattice orientation for both cases.

Secondary ion mass spectrometry was used to identify the isotopic composition of the single-crystalline graphite samples. After isotope enrichment using the high-pressure/high-temperature method54,55, our graphite crystals are composed of 99.98% 12C and 0.02% 13C, which are the world-class purified graphite samples.

To realize the reversal of heat flow in our Tesla valve, we separately designed the forward and reverse structures by rotating one by 180° to the other. We fabricated mirror-symmetry geometry of both structures suitable for the symmetric heat dissipation in the pump–probe measurement. This results in resistances at the source/sink for the forward/reverse directions being not strictly the same. However, in designing our samples, we considered this difference and found its impact to be negligible.

First, we conducted both Monte Carlo and FEM simulations in simplified geometry, without the arm of the Tesla valve (that is, the bent channel). These simulations aim to examine the impact of geometry difference at the source and the sink. Extended Data Fig. 3 shows that phonon travel times (in Monte Carlo) are the same in both samples, as well as the decay times (in FEM). In other words, in both ballistic and diffusive views of phonon transport, slight differences in geometry at the source and the sink do not make an observable difference.

Furthermore, the impact of contact sizes is smoothed out by the non-mathematically perfect heat source and sink areas. Indeed, in a theoretical system with a constant heat sink temperature at the ends of the bridge, contact areas would play an important role. But in the real experimental system, the heat source area is slowly heated up with a Gaussian laser beam to saturation, whereas the heat sink is a large conductive area that is far from a perfect isothermal boundary. Thus, the specific shape of the contact is less important. This is easy to see in the following FEM simulation, in which we changed the width of a graphite connector to the heat sink, as shown in Extended Data Fig. 8. Although the constriction should have had a stronger impact, we can notice a difference in the decay curves only if we sweep the width parameter (W) in a much wider range than in our experimental samples.

We used a non-contact, time-domain thermoreflectance method to analyse the thermal transport properties of our samples. A laser at 642 nm was pulsed with a duration of 2 μs to periodically pump heat through the Au transducer to the Tesla valve structures. A continuous-wave laser at 514 nm was used to investigate the change in reflectance (ΔR) on the metal transducer during heat dissipation. The variation of reflectance (ΔR) on the Au transducer was detected using the continuous-wave (probe) laser, which is equivalent to the variation of temperature by means of the coefficient of thermoreflectance (Cth), that is, ΔR/R = CthΔT. The TDTR signal thus shows the intrinsic decay of the heat injected by the pulsed (pump) laser through the Tesla valve structures in the microsecond range, as shown in Fig. 1c. For each measurement, the LabVIEW program fitted the decay curves with an exponential relation, ~exp(−t/τ), and collected the last 200 thermal decay times (τ). The corresponding τ in each measurement was taken from the mean of its Gaussian distribution, as shown in Fig. 1d. Before the measurement at each temperature, we spent sufficient time for the temperature stabilization with an uncertainty of <0.1 K to ensure thermal equilibrium between the sample and the holder in the cryostat. Moreover, the decay time and the thermal conductivity were experimentally verified to be independent of the probe laser power heating, even with a more than doubled laser intensity (Extended Data Fig. 9). We assume that the thermal conductivity is also reasonably expected to be independent of the pump laser power, which only induces temporal heat pulse and is more experimentally difficult to be precisely characterized. At each given temperature, we moved the lasers back and forth on the forward and reverse structures and performed repeated measurements up to ten times for each structure (forward or reverse) and the error bars shown in the plots were calculated from their standard deviations.

In the μ-TDTR measurement, we first quantify the heat dissipation through the sample by means of the thermal decay time, τ. The thermal conductivity of the Tesla valve is then extracted using FEM simulation implemented in COMSOL Multiphysics.

To reproduce the experiment, we built 3D models with the same geometries as the measured Tesla valve in forward and reverse directions, as shown in Extended Data Fig. 10. The pump laser centred on the gold transducer in the experiment is simulated by a Gaussian pulse heat source with the same pulse duration of 2 μs. We set the boundary condition of the heat sinks at the cryostat temperature. And the boundary condition of the edges of the Tesla valve structure is considered adiabatic, at which the thermal radiation loss of the suspended structure was estimated by the Stefan–Boltzmann law and concluded to be neglected, as well as the convection effect56,57.

Heat conduction through the entire model follows the heat diffusion equation:

in which ρ is the density, Cp is the heat capacity at constant pressure and κ is the thermal conductivity. These parameters are taken from the corresponding material values for the temperature of interest58,59. Thermal boundary conductance between the metal transducer and graphite has a negligible effect on our simulation results, as was investigated in our previous work10. We also consider the anisotropic thermal property in graphite by including the out-of-plane thermal conductivity from the literature59. Thus, the in-plane thermal conductivity of the Tesla valve structure is the only fitting parameter in the simulation. For instance, by sweeping a series of thermal conductivities of the graphite Tesla valve in the forward direction at 45 K, we obtain their corresponding decay curves, as plotted in Fig. 1c. Therefore, the corresponding thermal conductivity value can be obtained through the decay curve that fits our experimental data the best.

At 300 K, we obtained the thermal conductivity for the graphite Tesla valve comparable with the literature values. Our values are between those of the bulk graphite59 and the 1.2-μm-wide graphite ribbon with the same isotope contents60, as shown in Extended Data Fig. 6b. As the temperature decreases, these values are much lower than that of the bulk owing to the stronger size effects from the microscale geometry.

Our Monte Carlo algorithm traces phonon wave packets through the three-dimensional model of the structure that is assumed to be at temperature T. The phonons are emitted at the hot side and absorbed at the cold side. The phonon frequency (ω) and group velocities are obtained from the Planck distribution at the given temperature and phonon dispersion61. The dispersion was assumed to be isotropic and consisted of three acoustic branches. The temperature variations in the structure are assumed to be too small to cause substantial shift of the Planck distribution.

The Umklapp scattering events are simulated as scattering in a random direction that occurs once the time since the previous diffuse scattering event exceeds time t = −ln(r)τR, in which r is a random value between zero and one and τR is the relaxation time. The relaxation time is modelled as:

in which B = 3.18 × 10−25 sK−3 and Tdeb = 1,000 K (ref. 62). The impurity scattering is assumed negligible owing to the high purity of our sample. Normal-scattering events are not modelled owing to the nature of the current phonon Monte Carlo algorithm in the single-mode relaxation-time approximation. Thus, the simulations include no hydrodynamic effects and aim to demonstrate the parity of forward and reverse structures in the absence of phonon hydrodynamics. Modelling of boundary scattering depends on the surface roughness (σ), phonon wavelength (λ) and the incident angle to the surface (α). The probability (p) of boundary scattering to be specular or diffuse is evaluated for each scattering according to Soffer’s equation:

Then, a random value (r) is drawn to decide whether the scattering event is specular (r < p) or diffuse (r > p). In the case of diffuse scattering, the scattering-angle distribution follows the Lambert cosine law. The surface roughness was assumed to be 2 nm for side walls and 0.2 nm for top and bottom walls. Although simulations are performed for thousands of phonons, the maps of phonon trajectories show the paths of only a few first phonons, for clarity. The heat-flux maps show the in-plane flux defined as \(Q={({Q}_{x}^{2}+{Q}_{y}^{2})}^{0.5}\), in which Qx and Qy are the sums of the qx and qy components summed for all phonons integrated over time. The thermal maps are calculated as a sum of the energy (ħω) of all phonons at a given pixel integrated over time. The thermal conductivity is calculated as:

in which kb is the Boltzmann constant, ω(q) is the phonon frequency and v(q) is the group velocity at the wavevector q. The phonon relaxation time τR is measured using the Monte Carlo simulation as the average time between diffuse scattering events.

A previous experimental report has shown that the thermal conductivity of free-standing graphene exceeds 2,500 Wm−1 K−1 at room temperature63. However, when supported by copper, the thermal conductivity experiences a substantial reduction, decreasing to 370 Wm−1 K−1. The presence of a supporting substrate induces notable phonon–substrate interactions, thereby influencing the study of intrinsic thermal properties in the sample. Consequently, to circumvent substantial heat dissipation to the supporting substrate and ensure that heat is exclusively dissipated through the structures of interest, the fabrication of suspended graphite Tesla valves became imperative for the precise investigation of thermal rectification in this study.

However, microscale suspended structures unavoidably experience mechanical instability owing to residual stresses inherent in the fabrication process. At increased temperatures, the mechanical instability introduces an unavoidable discrepancy between successive measurements, resulting in discrepancies in the obtained thermal conductivity values at higher temperatures, as illustrated by the error bars in Fig. 4. It is a common challenge acknowledged in high-temperature measurements of suspended systems64, particularly for materials characterized by weak van der Waals interaction65, such as graphite in this study.

Despite these challenges, it is noteworthy that the diodicity values exhibit a consistent variation along the baseline of 1 from 70 to 300 K. To further validate our observations, further experiments were conducted on graphite Tesla valve samples S2 and S3. As depicted in Fig. 4, the findings from samples S2 and S3 confirm the high-temperature trends observed in sample S1, as discussed in the main text. By contrast, within the hydrodynamic temperature range of 25–60 K, at which the occurrence of thermal rectification is anticipated, we observed diodicity values consistently surpassing 1, extending beyond the margin of error.

All the data in the experiments and analysis that support the findings of this study are included in the main paper and Methods. Source data are provided with this paper.

The Monte Carlo FreePATHS algorithm is available at https://github.com/anufrievroman/freepaths.git (ref. 46).

Tesla, N. Valvular conduit. US patent 1,329,559 (1920).

de Vries, S. F., Florea, D., Homburg, F. & Frijns, A. Design and operation of a Tesla-type valve for pulsating heat pipes. Int. J. Heat Mass Transf. 105, 1–11 (2017).

Article ADS Google Scholar

Jin, Z.-j, Gao, Z.-x, Chen, M.-r & Qian, J.-y Parametric study on Tesla valve with reverse flow for hydrogen decompression. Int. J. Hydrogen Energy 43, 8888–8896 (2018).

Article ADS CAS Google Scholar

Porwal, P. R., Thompson, S. M., Walters, D. K. & Jamal, T. Heat transfer and fluid flow characteristics in multistaged Tesla valves. Numer. Heat Transf. A Appl. 73, 347–365 (2018).

Article ADS CAS Google Scholar

Nguyen, Q. M., Abouezzi, J. & Ristroph, L. Early turbulence and pulsatile flows enhance diodicity of Tesla’s macrofluidic valve. Nat. Commun. 12, 2884 (2021).

Article ADS CAS Google Scholar

Purwidyantri, A. & Prabowo, B. A. Tesla valve microfluidics: the rise of forgotten technology. Chemosensors 11, 256 (2023).

Article CAS Google Scholar

Lee, S., Broido, D., Esfarjani, K. & Chen, G. Hydrodynamic phonon transport in suspended graphene. Nat. Commun. 6, 6290 (2015).

Article ADS CAS Google Scholar

Huberman, S. et al. Observation of second sound in graphite at temperatures above 100 K. Science 364, 375–379 (2019).

Article ADS CAS Google Scholar

Jeong, J., Li, X., Lee, S., Shi, L. & Wang, Y. Transient hydrodynamic lattice cooling by picosecond laser irradiation of graphite. Phys. Rev. Lett. 127, 085901 (2021).

Article ADS CAS Google Scholar

Huang, X. et al. Observation of phonon Poiseuille flow in isotopically purified graphite ribbons. Nat. Commun. 14, 2044 (2023).

Article ADS CAS Google Scholar

Krishna Kumar, R. et al. Superballistic flow of viscous electron fluid through graphene constrictions. Nat. Phys. 13, 1182–1185 (2017).

Article CAS Google Scholar

Sulpizio, J. A. et al. Visualizing Poiseuille flow of hydrodynamic electrons. Nature 576, 75–79 (2019).

Article ADS CAS Google Scholar

Ku, M. J. et al. Imaging viscous flow of the Dirac fluid in graphene. Nature 583, 537–541 (2020).

Article ADS CAS Google Scholar

Cepellotti, A. et al. Phonon hydrodynamics in two-dimensional materials. Nat. Commun. 6, 6400 (2015).

Article ADS CAS Google Scholar

Guo, Y. & Wang, M. Phonon hydrodynamics and its applications in nanoscale heat transport. Phys. Rep. 595, 1–44 (2015).

Article ADS MathSciNet Google Scholar

Beardo, A. et al. Observation of second sound in a rapidly varying temperature field in Ge. Sci. Adv. 7, eabg4677 (2021).

Article ADS CAS Google Scholar

Ding, Z. et al. Observation of second sound in graphite over 200 K. Nat. Commun. 13, 285 (2022).

Article ADS CAS Google Scholar

Machida, Y. et al. Observation of Poiseuille flow of phonons in black phosphorus. Sci. Adv. 4, eaat3374 (2018).

Article ADS Google Scholar

Martelli, V., Jiménez, J. L., Continentino, M., Baggio-Saitovitch, E. & Behnia, K. Thermal transport and phonon hydrodynamics in strontium titanate. Phys. Rev. Lett. 120, 125901 (2018).

Article ADS CAS Google Scholar

Shang, M.-Y., Zhang, C., Guo, Z. & Lü, J.-T. Heat vortex in hydrodynamic phonon transport of two-dimensional materials. Sci. Rep. 10, 8272 (2020).

Article ADS CAS Google Scholar

Guo, Y., Zhang, Z., Nomura, M., Volz, S. & Wang, M. Phonon vortex dynamics in graphene ribbon by solving Boltzmann transport equation with ab initio scattering rates. Int. J. Heat Mass Transf. 169, 120981 (2021).

Article CAS Google Scholar

Zhang, C., Chen, S. & Guo, Z. Heat vortices of ballistic and hydrodynamic phonon transport in two-dimensional materials. Int. J. Heat Mass Transf. 176, 121282 (2021).

Article Google Scholar

Guo, Y. & Wang, M. Heat transport in two-dimensional materials by directly solving the phonon Boltzmann equation under Callaway’s dual relaxation model. Phys. Rev. B 96, 134312 (2017).

Article ADS Google Scholar

Ding, Z. et al. Phonon hydrodynamic heat conduction and Knudsen minimum in graphite. Nano Lett. 18, 638–649 (2018).

Article ADS MathSciNet CAS Google Scholar

Li, X. & Lee, S. Crossover of ballistic, hydrodynamic, and diffusive phonon transport in suspended graphene. Phys. Rev. B 99, 085202 (2019).

Article ADS CAS Google Scholar

Guo, Y. et al. Size effect on phonon hydrodynamics in graphite microstructures and nanostructures. Phys. Rev. B 104, 075450 (2021).

Article ADS CAS Google Scholar

Chen, G. Non-Fourier phonon heat conduction at the microscale and nanoscale. Nat. Rev. Phys. 3, 555–569 (2021).

Article CAS Google Scholar

Liao, B. (ed.) Nanoscale Energy Transport 2053–2563 (IOP, 2020).

Ghosh, K., Kusiak, A. & Battaglia, J.-L. Phonon hydrodynamics in crystalline materials. J. Phys. Condens. Matter 34, 323001 (2022).

Article CAS Google Scholar

Li, B., Wang, L. & Casati, G. Thermal diode: rectification of heat flux. Phys. Rev. Lett. 93, 184301 (2004).

Article ADS Google Scholar

Chang, C. W., Okawa, D., Majumdar, A. & Zettl, A. Solid-state thermal rectifier. Science 314, 1121–1124 (2006).

Article ADS CAS Google Scholar

Martínez-Pérez, M. J., Fornieri, A. & Giazotto, F. Rectification of electronic heat current by a hybrid thermal diode. Nat. Nanotechnol. 10, 303–307 (2015).

Article ADS Google Scholar

Shrestha, R. et al. Dual-mode solid-state thermal rectification. Nat. Commun. 11, 4346 (2020).

Article ADS CAS Google Scholar

Zhang, Y. et al. Simultaneous electrical and thermal rectification in a monolayer lateral heterojunction. Science 378, 169–175 (2022).

Article ADS CAS Google Scholar

Malik, F. K. & Fobelets, K. A review of thermal rectification in solid-state devices. J. Semicond. 43, 103101 (2022).

Article ADS CAS Google Scholar

Wang, H. et al. Experimental study of thermal rectification in suspended monolayer graphene. Nat. Commun. 8, 15843 (2017).

Article ADS CAS Google Scholar

Kasprzak, M. et al. High-temperature silicon thermal diode and switch. Nano Energy 78, 105261 (2020).

Article CAS Google Scholar

Desmarchelier, P., Tanguy, A. & Termentzidis, K. Thermal rectification in asymmetric two-phase nanowires. Phys. Rev. B 103, 014202 (2021).

Article ADS CAS Google Scholar

Wirtz, L. & Rubio, A. The phonon dispersion of graphite revisited. Solid State Commun. 131, 141–152 (2004).

Article ADS CAS Google Scholar

Lindsay, L., Broido, D. & Mingo, N. Flexural phonons and thermal transport in multilayer graphene and graphite. Phys. Rev. B 83, 235428 (2011).

Article ADS Google Scholar

Schelling, P. & Keblinski, P. Thermal expansion of carbon structures. Phys. Rev. B 68, 035425 (2003).

Article ADS Google Scholar

Lee, S., Li, X. & Guo, R. Thermal resistance by transition between collective and non-collective phonon flows in graphitic materials. Nanoscale Microscale Thermophys. Eng. 23, 247–258 (2019).

Article ADS CAS Google Scholar

Huang, X., Guo, Y., Volz, S. & Nomura, M. Mapping phonon hydrodynamic strength in micrometer-scale graphite structures. Appl. Phys. Express 15, 105001 (2022).

Article ADS CAS Google Scholar

Huang, X. et al. Coherent and incoherent impacts of nanopillars on the thermal conductivity in silicon nanomembranes. ACS Appl. Mater. Interfaces 12, 25478–25483 (2020).

Article CAS Google Scholar

Anufriev, R. & Nomura, M. Ray phononics: thermal guides, emitters, filters, and shields powered by ballistic phonon transport. Mater. Today Phys. 15, 100272 (2020).

Article Google Scholar

FreePATHS - free phonon and thermal simulator. GitHub https://github.com/anufrievroman/freepaths (2024).

Ravichandran, N. K. & Broido, D. Phonon-phonon interactions in strongly bonded solids: selection rules and higher-order processes. Phys. Rev. X 10, 021063 (2020).

CAS Google Scholar

Klarbring, J., Hellman, O., Abrikosov, I. A. & Simak, S. I. Anharmonicity and ultralow thermal conductivity in lead-free halide double perovskites. Phys. Rev. Lett. 125, 045701 (2020).

Article ADS CAS Google Scholar

Guyer, R. & Krumhansl, J. Thermal conductivity, second sound, and phonon hydrodynamic phenomena in nonmetallic crystals. Phys. Rev. 148, 778 (1966).

Article ADS CAS Google Scholar

Kim, W. Strategies for engineering phonon transport in thermoelectrics. J. Mater. Chem. C 3, 10336–10348 (2015).

Article CAS Google Scholar

Maznev, A., Every, A. & Wright, O. Reciprocity in reflection and transmission: what is a ‘phonon diode’? Wave Motion 50, 776–784 (2013).

Article ADS Google Scholar

Geurs, J. et al. Rectification by hydrodynamic flow in an encapsulated graphene Tesla valve. Preprint at https://arxiv.org/abs/2008.04862 (2020).

Hu, J., Ruan, X. & Chen, Y. P. Thermal conductivity and thermal rectification in graphene nanoribbons: a molecular dynamics study. Nano Lett. 9, 2730–2735 (2009).

Article ADS CAS Google Scholar

Taniguchi, T. & Yamaoka, S. Spontaneous nucleation of cubic boron nitride single crystal by temperature gradient method under high pressure. J. Cryst. Growth 222, 549–557 (2001).

Article ADS CAS Google Scholar

Taniguchi, T. & Watanabe, K. Synthesis of high-purity boron nitride single crystals under high pressure by using Ba–BN solvent. J. Cryst. Growth 303, 525–529 (2007).

Article ADS CAS Google Scholar

Pope, A., Zawilski, B. & Tritt, T. Description of removable sample mount apparatus for rapid thermal conductivity measurements. Cryogenics 41, 725–731 (2001).

Article ADS CAS Google Scholar

Maire, J. Thermal Phonon Transport in Silicon Nanosturctures. PhD thesis, Univ. Tokyo (2015).

Nihira, T. & Iwata, T. Temperature dependence of lattice vibrations and analysis of the specific heat of graphite. Phys. Rev. B 68, 134305 (2003).

Article ADS Google Scholar

Ho, C. Y., Powell, R. W. & Liley, P. E. Thermal conductivity of the elements. J. Phys. Chem. Ref. Data 1, 279–421 (1972).

Article ADS CAS Google Scholar

Huang, X. et al. Super-ballistic width dependence of thermal conductivity in graphite nanoribbons and microribbons. Nanomaterials 13, 1854 (2023).

Article CAS Google Scholar

Monteverde, U. et al. Under pressure: control of strain, phonons and bandgap opening in rippled graphene. Carbon 91, 266–274 (2015).

Article CAS Google Scholar

Alofi, A. & Srivastava, G. Thermal conductivity of graphene and graphite. Phys. Rev. B 87, 115421 (2013).

Article ADS Google Scholar

Cai, W. et al. Thermal transport in suspended and supported monolayer graphene grown by chemical vapor deposition. Nano Lett. 10, 1645–1651 (2010).

Article ADS CAS Google Scholar

Aliane, A. et al. Mechanical modeling and characterization of suspended cooled silicon bolometers for sub-millimeter and millimeter waves polarization detection. Sens. Actuators A Phys. 296, 254–264 (2019).

Article ADS CAS Google Scholar

Wang, M. C. et al. Mechanical instability driven self-assembly and architecturing of 2D materials. 2D Mater. 4, 022002 (2017).

Article Google Scholar

Asheghi, M., Touzelbaev, M. N., Goodson, K. E., Leung, Y. K. & Wong, S. S. Temperature-dependent thermal conductivity of single-crystal silicon layers in SOI substrates. J. Heat Transf. 120, 30–36 (1998).

Article CAS Google Scholar

Tang, J. et al. Holey silicon as an efficient thermoelectric material. Nano Lett. 10, 4279–4283 (2010).

Article ADS CAS Google Scholar

Download references

This work was supported by the JSPS Grants-in-Aid for Scientific Research (grant no. 21H04635), JST SICORP EIG CONCERT-Japan (grant no. JPMJSC22C6) and ALCA-Next (grant no. JPMJAN23E3), Japan Science and Technology Agency. K.W. and T.T. acknowledge support from the JSPS KAKENHI (grant nos. 20H00354, 21H05233 and 23H02052) and World Premier International Research Center Initiative (WPI), MEXT, Japan. We thank X. Luo, Y. Kosevich, Y. Wu, X. Wu and L. Yang for discussions, R. Yanagisawa for the help with sample fabrication and P. Gassmann for the improvement of Monte Carlo simulations.

Institute of Industrial Science, The University of Tokyo, Tokyo, Japan

Xin Huang, Roman Anufriev, Laurent Jalabert, Sebastian Volz & Masahiro Nomura

Université de Lyon, INSA Lyon, CNRS, CETHIL, UMR5008, Villeurbanne, France

Roman Anufriev

LIMMS, CNRS-IIS IRL 2820, The University of Tokyo, Tokyo, Japan

Laurent Jalabert, Sebastian Volz & Masahiro Nomura

Research Center for Electronic and Optical Materials, National Institute for Materials Science, Tsukuba, Japan

Kenji Watanabe

Research Center for Materials Nanoarchitectonics, National Institute for Materials Science, Tsukuba, Japan

Takashi Taniguchi

School of Energy Science and Engineering, Harbin Institute of Technology, Harbin, People’s Republic of China

Yangyu Guo

School of Physical Science and Technology, Southwest Jiaotong University, Chengdu, People’s Republic of China

Yuxiang Ni

You can also search for this author in PubMed Google Scholar

You can also search for this author in PubMed Google Scholar

You can also search for this author in PubMed Google Scholar

You can also search for this author in PubMed Google Scholar

You can also search for this author in PubMed Google Scholar

You can also search for this author in PubMed Google Scholar

You can also search for this author in PubMed Google Scholar

You can also search for this author in PubMed Google Scholar

You can also search for this author in PubMed Google Scholar

X.H. led the research with M.N., designed and fabricated the samples, performed the TDTR measurements, analysed the results and wrote the manuscript. R.A. performed Monte Carlo simulations and contributed to the sample design. R.A. and L.J. contributed to the TDTR experiments. K.W. and T.T. synthesized the isotopically purified graphite crystals. Y.G., Y.N. and S.V. contributed to interpreting the results. M.N. conceived the original idea, conceptualized the study, designed and supervised the overall research project, secured the necessary financing and was responsible for project administration. All authors contributed to discussing the results and manuscript revision.

Correspondence to Xin Huang or Masahiro Nomura.

The authors declare no competing interests.

Nature thanks the anonymous reviewers for their contribution to the peer review of this work. Peer reviewer reports are available.

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

a, Scanning electron microscopy image of the graphite sample flake used for EBSD. b, Scheme of EBSD measurement. The EBSD pattern carries information about the crystal orientation in a sample as indicated by the inverse pole figures for the normal direction (c), the transverse direction (d) and the reference direction (e).

Suspended graphite Tesla valves vertical (a) and parallel (b) to the heat-dissipation directions.

Monte Carlo (a) and FEM (b) simulations in the main channels of the Tesla valve with geometry difference at the sources and the sinks.

Source Data

Thermal energy distribution in the Tesla valve in both the forward and reverse directions at 40 K using Monte Carlo simulations.

Temperature-dependent thermal conductivity and diodicity of graphite Tesla valve samples S1 (a) and S2 (b) in the forward (cyan triangles) and reverse (red triangles) directions. Thermal conductivity (top) and diodicity (bottom) of samples S1 (c) and S2 (d) from 20 to 60 K in linear scale.

Source Data

a, Temperature-dependent thermal conductivity of the silicon Tesla valve in the forward (open cyan triangles) and reverse (open red triangles) directions with a thickness of 145 nm (this work) and silicon membranes with thicknesses of 420 nm (ref. 66) and 100 nm (ref. 67). b, Temperature-dependent thermal conductivity of graphite Tesla valve sample S3 in both forward and reverse directions (this work), bulk pyrolytic graphite59 and 1.2-μm-wide graphite ribbon60.

Source Data

a, Thermal conductivity and diodicity as a function of temperature from 20 to 60 K (hydrodynamic temperature range of graphite) in normal scale. Inset, optical microscope image of Tesla valves in the forward (top) and reverse (bottom) directions patterned on a 145-nm-thick silicon membrane. The lengths and widths of the channels are the same as those of the graphite Tesla valve. b, Thermal conductivity and diodicity of a supplementary silicon sample.

Source Data

FEM simulations in graphite ribbon with different contact areas at the heat sinks.

Source Data

Measured thermal conductivity of the graphite Tesla valve in the forward and reverse directions as a function of probe laser power in the μ-TDTR system at 300 K.

Source Data

Illustration of heat dissipation in the graphite Tesla valve in both the forward (left) and reverse (right) directions in FEM models.

Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.

Reprints and permissions

Huang, X., Anufriev, R., Jalabert, L. et al. A graphite thermal Tesla valve driven by hydrodynamic phonon transport . Nature (2024). https://doi.org/10.1038/s41586-024-08052-1

Download citation

Received: 02 March 2024

Accepted: 16 September 2024

Published: 16 October 2024

DOI: https://doi.org/10.1038/s41586-024-08052-1

Anyone you share the following link with will be able to read this content:

Sorry, a shareable link is not currently available for this article.

Provided by the Springer Nature SharedIt content-sharing initiative