A bottom-up route to enhance thermoelectric figures of merit in graphene nanoribbons.

We propose a hybrid nano-structuring scheme for tailoring thermal and thermoelectric transport properties of graphene nanoribbons. Geometrical structuring and isotope cluster engineering are the elements that constitute the proposed scheme. Using first-principles based force constants and Hamiltonians, we show that the thermal conductance of graphene nanoribbons can be reduced by 98.8% at room temperature and the thermoelectric figure of merit, ZT, can be as high as 3.25 at T = 800 K. The proposed scheme relies on a recently developed bottom-up fabrication method, which is proven to be feasible for synthesizing graphene nanoribbons with an atomic precision.

Despite the fact that it has been more than 150 years since the discovery of thermoelectric (TE) effect, its applications are still limited1. This limitation stems from the apparently incompatible demands for high TE efficiency which is related to the material specific dimensionless figure of merit defined as ZT = S2GT/κ with S being the Seebeck coefficient, G being the electrical conductance, T being temperature and κ = κel + κph is the thermal conductance consisting of electron and phonon contributions. A good TE material needs to have a high electrical conductivity like in a metal, a high Seebeck coefficient like in an insulator and a low thermal conductivity like in a glass2. Meeting these diverse, if not contradictory, attributes in a single system is an extremely difficult challenge. Nano-scale fabrication techniques, on the other hand, enable the tailoring of the material properties to a previously unprecedented extend. By reducing the cross section of a nanowire, its TE performance can be enhanced significantly due to the confinement of electrons3. Nano-structuring also offers ways to control phonon transport by tuning phonon dispersions and introducing scatterers specially designed to scatter the dominant phonon wavelengths, thus enhance ZT45678910111213141516. Considering the fact that low frequency acoustic phonons dominate phonon conduction, the size of the phonon scatterers plays an important role in lowering thermal conductivity17181920. It is worth mentioning that the knowledge of ballistic transport properties is a necessary starting point for investigating TE efficiency2122. On the other hand, addressing the scattering mechanisms is not only required for an accurate prediction, but engineering of scatterers is a promising way to enhance the TE performance23. Indeed, state-of-the-art TE materials design relies on engineering of scattering mechanisms for heat and charge carriers so as to approach the phonon glass-electron crystal limit24.

Graphene possesses a wide range of superlative material properties25, including the record value for thermal conductivity262728293031323334353637. Still, it is possible to tailor the thermal transport properties of graphene and graphene nanoribbons' (GNRs) by nano-structuring techniques such as edge roughness1113, defect engineering1338, isotope engineering2039 and introducing periodic nano-holes40. Bottom-up fabrication of GNRs using surface-assisted coupling of molecular precursors is not only a promising technique to obtain atomically precise GNRs with predefined geometries and sub-nanometer widths41, but it also enables a hybrid scheme for TE efficiency by combining geometrical structuring with isotope clustering. This combination enables a reduction of the thermal conductivity by up to two orders of magnitude. Also, electronic mini-bands are formed due to the chevron geometry, approaching the Mahan-Sofo condition42. The electronic bandwidths are compatible with the optimal values defined by Zhou et al23. As a result, the Seebeck coefficient and the power factor are substantially enhanced and thus ZT values up to 3.25 are possible.

We consider two types of GNRs which have been synthesized by Cai et al41. The first type has an armchair edge shape with 7 dimers per unit cell and a straight geometry, s-GNR, that was obtained from 10,10′-dibromo-9-9′-bianthryl precursor monomers. The second GNR type is obtained from tetraphenyl-triphenylene monomers and has a chevron-type geometry, c-GNR. Both GNRs have well-defined geometries and their edges are passivated with hydrogen. (see Figure 1a) The width of a ribbon, w, is defined as the total area of hexagons in its unit cell, divided by the unit cell length (w = 0.62 nm and 0.93 nm for s-GNR and c-GNR, respectively), and we use the interlayer distance of graphite as the ribbon height. We employ the density functional tight binding (DFTB) approach to calculate the electronic and vibrational properties of the GNRs43. Having obtained the force constant matrices and the electronic Hamiltonians from DFTB simulations, we use atomistic Green's functions (AGF)44 and nonequilibrium Green's functions (NEGF)45 for calculating the phonon and electron transport properties, respectively. (see Supplementary Information for details)

Figure 1

Structural aspects, isotope distribution types, and ballistic transmission spectra of bottom-up fabricated graphene nanoribbons (GNRs).

(a) Structures of precursors and corresponding straight and chevron type graphene nanoribbons, s-GNR and c-GNR. (b) Heavy isotopes can be distributed at the atomic or precursor level. Grey, 12C; black, 14C; white H. Ballistic phonon (c) and electron (d) transmission spectra are plotted for both GNR types, where mini band formation due to the geometry of c-GNR is evident.

Results

Ballistic phonon and electron transmission spectra of s-GNR and c-GNR are plotted in Figures 1c and 1d, respectively. The first observation is that the phonon dispersions of the s-GNR and c-GNR are significantly different due to their geometries. The chevron shape gives rise to the formation of mini-bands with small group velocities and numerous phonon energy gaps arise (see Figure 1c and the Supplementary Information). Acoustic phonons, which have large dispersions and very long mean free paths in straight GNRs, have significantly smaller group velocities in c-GNRs. As a result, ballistic phonon transmission is significantly altered (Figure 1c) and ballistic lattice thermal conductance per cross section area, κph/A, is reduced by 69% at room temperature (see Figure 2 and Table I). More importantly, the phonon density of states (DOS) of c-GNR posseses many van Hove singularities in the entire spectrum giving rise to enhanced scattering rates upon the inclusion of scatterers.

Figure 2

Phonon transport through isotopically disordered GNRs.

Thermal conductance per cross section area (κ/A) are plotted as functions of temperature in (a) and ribbon length (L) in (b), for s-GNR (blue) and c-GNR (red). Thermal conductance is minimum for maximal disorder (d = 50%) for a given distribution type. Precursor distribution always yields lower κ/A for a given isotope density for s-GNR. For c-GNR, precursor distribution yields larger κ/A at low density and short L. For long systems, even a low density of heavy precursors give rise to stronger suppression of phonon transport than that of atomic distribution of isotopes. Solid lines in (b) depict κ/A in the absence of 14C isotopes.

Table 1

Suppression of lattice thermal transport. κ/A for pristine and 5 μm long GNRs with 50% 14C isotopes with atomic and precursor distribution of isotopes. T = 300 K, 500 K and 800 K. κ/A is given in units of nW K−1 nm−2

		pristine	atomic (50%)	precursor (50%)
s-GNR	300 K	5.002	0.566	0.460
	500 K	7.187	0.585	0.477
	800 K	8.935	0.592	0.484
c-GNR	300 K	1.571	0.131	0.060
	500 K	2.310	0.131	0.060
	800 K	2.914	0.132	0.060

We introduce 14C isotopes in order to suppress phonon transport without damaging the electronic quality with (i) random atomic distribution of isotopes, and (ii) random distribution of heavy precursors (Firgure 1b). In the first case, each 14C atom is a scattering center, while in the second case the precursors consist of only 12C or only 14C. Clusters of impurities are introduced in order to overcome the alloy limit. This is due to the fact that scattering rates increase with the size of the scatterer for long wavelength phonons, which dominate lattice thermal conduction. For both s-GNR and c-GNR, we consider both atomic and precursor distributions with 14C densities being d = 10% and 50%.

First, we analyze the effects of different isotope distributions on the phonon transport properties of the straight s-GNR (Figure 2a–b). Using the AGF technique, we calculate transmission spectra at various lengths for ensembles consisting of at least 25 samples. The phonon mean-free-paths (MFPs) are obtained using . Here, L is the ribbon length, stands for ensemble averaged transmission value, is the ballistic transmission. For the sake of comparison, we also estimate the phonon MFPs from transmission spectra of individual scatterers, by making use of scaling theory as Here, Nuc is the number of atoms in the unit cell, luc is the length of the unit cell with i being the atomic index, and is the transmission spectrum of GNR when only the ith atom is isotopically different in the entire system. The scaling approximation works well when d is small so that interference between different scattering events can be neglected. agrees well with for d ≤ 10%. For d > 10%, scaling theory often predicts longer mean-free-paths (MFPs), because it disregards multiple scattering effects. That is, the relation ∝ 1/d does not hold for d > 10%. Comparing the values of for different isotope densities, it is found that d = 10% always yields longer than d = 50%, when the distribution type (atomic or precursor is kept the same), which is also predicted from Equation 1. On the other hand, different distribution types, when d is kept constant, give rise to opposite behavior at low and high frequency regions of the spectrum. At low frequency, the atomic distribution generally yields longer than the precursor distribution, while the opposite holds at high frequencies1920. (see Supplementary Information) In the Rayleigh limit, the scattering cross section increases with , where N is the number of atoms in the scattering center. At high frequencies, on the other hand, the scattering cross section scales as 46. One route to suppress lattice thermal conduction is the so called nanoparticle-in-alloy approach, which incorporates clusters of scatterers together with alloying and aims to shorten mean-free-paths in the entire phonon spectrum17. In our approach, effective suppression of high-frequency phonons is realized by combining the effects of geometrical structuring, i.e chevron geometry with clustered isotopes. Similar to s-GNR, low frequency phonons in c-GNR are scattered more effectively by precursor distribution of isotopes while atomic scatterers suppress high energy phonons more strongly (see Supplementary Information). Combined with the effect of geometry on phonon dispersions, is always shorter than 100 nm for ω > 50 cm−1 at 50% isotope density with precursor distribution, and when the c-GNR is longer than 1 μm these phonons have negligibly small contribution to heat transport.

Phonon thermal conductance is calculated from the ensemble averaged transmission values as where f is the Bose-Einstein distribution function. In Figure 2, the phonon thermal conductance per cross section area (κph/A) is plotted as a function of temperature and GNR length (L), for pristine and isotopically engineered GNRs having atomic and precursor isotope distributions with d = 10% and 50% for s-GNR (blue) and c-GNR (red). For s-GNR, d = 50% always yields lower κph than d = 10% for both distribution types and atomic distribution gives rise to higher κph. For c-GNR, on the other hand, the atomic distribution of isotopes yields higher κph for L > 0.08 μm independent of d, and higher d results in lower κph for a given distribution type. The room temperature behavior of κph for c-GNRs with a precursor distribution at different L with d = 10% is a consequence of the frequency and length dependence of phonon transmission on distribution type (Figure 2b). κph is highest with dprec = 10% at L < 0.01 μm, where datom = 10% yields the highest κph for L > 0.02 μm. For L > 0.08 μm, dprec = 10% results in κph lower than datom = 50%. This is because (i) at short L, phonons with ω > 50 cm−1 have appreciable contribution to κph and they have longer for precursor distribution; (ii) for long L, phonons with ω > 50 cm−1 are strongly suppressed independent of the distribution type for both d (i.e. only low frequency phonons contribute to κph) and is shorter for precursor distribution of isotopes. As a result, the precursor distribution of isotopes is highly efficient in reducing the thermal conductivity of c-GNR.

Comparing the area normalized phonon thermal conductance values in Table I, we find that the ballistic value of κph/A is reduced by 69% due to chevron geometry. Isotopic disorder induces suppressions of 97% (89%) and 98.8% (91%) for c-GNR (s-GNR) with precursor and atomic distributions, respectively. In other words, compared to the atomic distribution, κ/A is reduced by 54% (18%) for the precursor distribution in c-GNRs (s-GNRs). As a comparison, thermal conductivity of isotopically disordered two-dimensional graphene is lower than the isotopically pure sample by about a factor of 2 and in GNRs reductions by a factor of 3 with respect to the pure case were reported3947.

The geometrical aspects of GNRs are not only crucial for the suppression of phonon transport, but they also enable electronic band engineering for tailoring the TE efficiency. The formation of mini-bands in c-GNR narrows the dispersion of electrons participating in transport. It was shown by Mahan and Sofo that, if one disregards scatterings, a Diracdelta shaped resonance in the electronic DOS close to the Fermi energy constitutes the optimal electronic structure for TE performance42. That is, the sharper the resonance, the better the TE performance that can be. In c-GNR, the band widths range between 0.1 to 0.2 eV and an extremely high TE efficiency is predicted within the ballistic electron assumption (see Supplementary Information). On the other hand, for a realistic calculation it is necessary to address the scattering effects in the electronic transport as well. In this case, TE performance is largely determined by the dimensionality of the system23. Anderson disorder48, for example, resembles the scattering model where the relaxation time is inversely proportional to the DOS, and it was previously utilized by White et al. to explore the electron transport length-scales in carbon nanotubes49. In the simplest case of a quasi-one-dimensional wire with a single band having dispersion E(k) = –W/2 cos(ka), the electronic MFP is (E) = aW2/4σ2(1–(2E/W)2), where W is the band with, a is the unit cell length and σ is the standard deviation of onsite disorder. As pointed out by Zhou et al.23, when W approaches zero, the mean free path and the transmission vanish and ZT becomes zero. Therefore W > 2.4 kBT is required for optimum TE performance23.

Returning to the c-GNR, we perform NEGF transport calculations with Anderson disorder (details explained in the Supplementary Information). The random onsite energies are chosen from an energy range whose width, w, is set to be proportional to the temperature as . Having calculated the transmission spectra for each sample, we average over an ensemble of 25 samples to get and obtain the TE coefficients as a function of temperature (T), chemical potential (μ) and sample length (L) (see Supplementary Information). The maxima of electric conductance, G, is reduced at elevated temperatures due to enhanced electron scattering, which reduces the Seebeck coefficient, S, as well. On the other hand, the peak values in the power factor (P = S2G) close to the edges of the valence and conduction bands are not affected by temperature significantly. This is mainly due to the enhanced values of G with T, when μ is inside the band gap. This effect is specific to the valence and conduction bands, while P is suppressed at elevated temperatures for the rest of the bands. Around the charge neutrality point κ is dominated by κph, which is suppressed by scattering from isotopic precursors. Despite the fact that ZTmax is considerably reduced compared to the ballistic electron assumption, where ZT can be as high as 7 for long GNRs, ZT ≥ 2 is achieved at room temperature for 1 μm > L > 0.1 μm. Higher ZT is possible at higher temperatures for shorter c-GNRs, e.g. ZT = 3.25 is predicted at T = 800 K and  μm (Figures 3b and 4).

Figure 3

Thermoelectric figure of merit, ZT, for chevron-type GNR with randomly distributed heavy precursors.

ZT is plotted as a function of chemical potential (μ) and length (L) at T = 800 K, (a). The heavy precursor density is d = 50%. Local maxima appear close to the band edges, and the maximum value is obtained inside the band gap for all lengths. ZT for optimum system lengths are plotted at T = 300 K (blue), 500 K (green)and 800 K (red), (b). Electron density of states is depicted in gray. The system lengths are L = 430 nm, 140 nm and 70 nm, respectively. ZT = 3.25 is realized at 800 K.

Figure 4

Maximum ZT as a function of length at different temperatures.

Maximum ZT achievable when Anderson type disorder is introduced in the electronic Hamiltonian. The variation of onsite energies are set equal to the temperature, σ = kBT. Solid lines represent hole-like transport while the dashed lines are for electron-like charge carriers. The ZTmax values shown are realized inside the band gap, i.e. |μ| < 0.75 eV with respect to the mid-gap, except for the dotted curve (T = 300 K) when |μ| > 1.

Discussion

The length dependence of maximum ZT at different temperatures, ZTmax, displays the trade off due to electron and phonon scatterings in TE transport (Figure 4). In contrast to the pristine electron prediction that ZTmax increases monotonically with length, electron scattering results in a substantial reduction with L. At T = 300 K, a high ZT is assured for a wide range of c-GNR lengths from 0.01 to 5 μm, while at T = 500 K and 800 K the ranges are narrower and the peak positions shift to shorter L. Stronger suppression of electronic transport at higher temperatures is the reason for the narrowing and the shift of the peak positions in the ZTmax curves. On the other hand, κph is almost insensitive to temperature for L ≥ 0.1 μm, which is due to extremely short phonon MFPs as a result of hybrid nanostructuring. P being robust to changes in temperature for both valence and conduction bands, and κph being saturated at relatively low temperatures are key aspects in the enhanced ZT values. As a final comment, s-GNRs with the same type of distribution of isotopes and Anderson disorder have ZT values smaller than 0.5 in the given range of L.

In summary, the combination of geometrical structuring and isotope engineering at the precursor level makes it possible to optimize both electronic and phononic transport properties of the considered carbon nanomaterials which have already been realized with existing fabrication technology. Noting that ZT > 3 is the goal for efficient thermoelectrics5051, isotopically engineered c-GNRs are good candidates for technological applications, especially at elevated temperatures. The proposed hybrid nanostructuring scheme is promising for efficient TE energy conversion and thermal management of nano-devices. This scheme is not limited to carbon based systems but it is applicable to low-dimensional structures in general.

Methods

In this work, electronic Hamiltonians and overlap matrices as well as the interatomic force constants are obtained using the density functional tight binding (DFTB) method as implemented in the DFTB+ sofware package5253. The electronic transport problem is treated using the non-equilibrium Green's function (NEGF) formalism45, while atomistic Green's function (AGF) technique is implemented for phonons44. Both electronic and phononic Green's functions are calculated recursively54. Thermoelectric coefficients are obtained from microscopic transport calculations55. More details about implementations of the methods are explained in the Supplementary Information.

Author Contributions

H.S. conceived the research. C.S. did the DFTB simulations. H.S. did the NEGF and AGF calculations. T.C. and G.C. oversaw all research phases. Everyone contributed to the writing of the paper.