Electron and phonon transport in twisted graphene nanoribbons

TL;DR

This paper investigates how torsional deformations affect electrical, thermal, and thermoelectric transport in graphene nanoribbons using advanced theoretical models.

Contribution

It introduces a comprehensive theoretical framework combining nonequilibrium Green's functions with ab-initio and empirical methods to analyze twisted graphene nanoribbons.

Findings

Twisting alters electron eigenchannels and phonon transmission.

Different twist configurations impact thermoelectric efficiency.

Transport properties are sensitive to the type and angle of deformation.

Abstract

We theoretically study the electrical, thermal and thermoelectric transport properties of graphene nanoribbons under torsional deformations. The modelling follows a nonequilibrium Green's function approach in the ballistic transport regime, describing the electrical and phononic properties through \textit{ab-initio} density functional theory and empirical interatomic potentials, respectively. We consider two different types of deformations, a continuous twist of a given angle applied to the nanoribbon, and two consecutive twists applied in opposite angular directions. The numerical results are carefully analysed in terms of spatially-resolved electron eigenchannels, polarization-dependent phonon transmission and thermoelectric figure-of-merit.

Figures

Equations

$${\mathcal{T}}^e(E)=Tr\left[{\mathbf{\Gamma }}_L(E){\mathbf{G}}_C^r(E){\mathbf{\Gamma }}_R(E){\mathbf{G}}_C^a(E)\right],$$

$$G^e(\mu )=\frac{2e^2}{h}\int {\mathcal{T}}^e(E)\left(-\frac{\partial f_0(E,\mu )}{\partial E}\right)dE,$$

$$G^{ph}(T)=\frac{\hslash }{2\pi }\int \omega {\mathcal{T}}^{ph}(\omega )\left(\frac{\partial n_0(\omega ,T)}{\partial T}\right)d\omega ,$$

$$ZT=\frac{S^2{G}^{e}T}{G^{ph}+{\kappa }^e},$$

$$S(\mu )=\frac{1}{eT}\frac{L_1(\mu )}{L_0(\mu )},$$

$${\kappa }^e(\mu )=\frac{1}{T}\left\{L_2(\mu )-\frac{L_1(\mu {)}^2}{L_0(\mu )}\right\},$$

$$L_m(\mu )=\frac{2}{h}{\int }_{-\infty }^{\infty }{\mathcal{T}}^e(E)(E-\mu {)}^m\left(-\frac{\partial f_0(E,\mu )}{\partial E}\right)dE.$$

Full text

††: J. Phys. D: Appl. Phys.

rrurali@icmab.es,

• March 2017

Keywords: Graphene nanoribbon, electron transport, thermal transport, atomistic simulation.

1 Introduction

The field of research in two-dimensional (2D) materials has been enjoying extraordinary growth during the past decade. This activity was triggered by pioneering works on graphene [1, 2, 3], a 2D semimetallic allotrope of carbon that turned out to be an exceptionally fertile ground for advancing frontiers of condensed matter physics [4, 5, 6, 7]. Today the family of 2D materials has become densely populated with manifold specimens [8, 9, 10]: it comprises monolayer materials made of a single element (phosphorene, borophene, germanene and silicene) and others featuring different atoms alternating in the same layer (boron nitride, transition metal dichalcogenides (TMDCs) and Mxenes).

Graphene consists of a hexagonal monolayer network of $sp^{2}$-hybridized carbon atoms whose properties were expected to be outstanding, based on theoretical predictions. Unique ballistic transport properties, very long mean free path at room temperature  [11], distinctive integral and half-integral quantum hall effect[12, 3], the highest known thermal conductivity [13] and high electron mobility [14] are among the most intriguing features of graphene. In particular, thanks to its electron mobility significantly higher than that of the widely-used Si, graphene has been envisioned as a candidate material for new generations of nano-electronic devices [15, 16].

From the experimental viewpoint, recent developments in fabrication techniques have made it possible to grow very narrow [17, 18] and atomically precise [19] graphene nanoribbons (GNRs), further stimulating interest in this 2D material. Nevertheless, like in any other real material, structural defects do exist and can dramatically alter the properties of graphene. Several experimental studies have shown, for instance, the occurrence of intrinsic or extrinsic defects in graphene [20, 21, 22, 23, 24]. Others demonstrated how the surface of fabricated graphene nanoribbons was not perfectly flat due to the presence of ripples, i.e. nanoscopic roughening of the 2D plane [25, 26].

These lattice imperfections have a strong influence on the electronic, optical, thermal, and mechanical properties of the material. As a matter of fact, many of the characteristics of technologically important materials such as the electrical conductance of semiconductors are governed by defects [27]. For this reason, defects are often deliberately introduced, by irradiation or chemical methods. This is the basis for the shaping of the material properties via engineering of its defects, a possibility which determined the need for a careful analysis of the relationship between the nature and amount of graphene defects and the deriving physical properties.

Several studies have, indeed, been devoted to the effect of structural defects on graphene properties, with major attention given to point defects, multiple vacancies and substitutional impurities [28, 29, 30]. Furthermore, it has been proved that the occurrence of ripples in the surface constitutes an intrinsic feature of 2D materials and its link with the mechanical stability of the layer has been explained [31, 4, 25]. For the specific case of GNRs, researchers have investigated the changes in thermal and electron transport due to edge terminations, shape, width and roughness [32, 33, 34, 35, 36]. Recently, interest is growing also on the their mechanical and electronic properties under different constraints such as mechanical stress [37, 38, 39, 40] From an application point of view, their flexibility could allow their employment in stretchable electronics and in the creation of new-generation devices for non-linear energy harvesting [41, 42, 43]. In this regard, of particular interest is twisting, as this is a mechanical deformation unique from GNRs due to their exceptional flexibility. Twisted nanoribbons have been fabricated encapsulated inside graphene nanotubes [44, 45] or by means of chemical etching. [46] The effect of twisting on their structural, [47, 48] electronic [49] and transport [50, 51, 52, 53] properties has been recently explored mostly from the theory side. Apart from showing the expectable degradation of the electronic and thermal transport properties, calculations performed at different twist angles have predicted an electromechanical switching effect induced by twisting and occurring in both armchair [51, 52] and zig-zag GNRs. [50]

In this work we investigate the electron and phonon transport properties of armchair graphene nanoribbons by considering two alternative torsional deformation processes: 1) a fixed twist of $180^{\circ}$ is applied along the longitudinal axis while the nanoribbon length is varied; 2) a twist of variable angle $\phi$ is applied for half the length of the nanoribbon followed by a twist of $-\phi$ over the remaining half. We study the evolution of the electrical and thermal conductance and analyse the results in terms of twist-induced structural deformation, spatially-resolved electron eigenchannels, polarization-dependent phonon transport and thermoelectric properties.

2 Computational methods

2.1 Electronic structure and electron transport

The electronic structure is calculated from first-principles, using the implementation of density-functional theory (DFT) of the Siesta package [54]. We account for the core electrons with norm-conserving pseudopotentials and use the Generalized Gradient Approximation (GGA) for the exchange-correlation energy. We have optimized a light, but accurate single-$\zeta$ polarized basis set to expand the one-electron wavefunction that allowed us dealing efficiently with relatively large systems.

We study an armchair graphene nanoribbon (GNR) of width 1.2 nm and variable length. The edge dangling bonds are passivated with hydrogen atoms. We sample the Brillouin zone with the $\Gamma$ point only, because of the very large size of the computational cell along the only periodic direction, taken to be parallel to the $z$-axis. We use a $1\times 1\times 4$ k-point mesh for the bulk calculations of the 3-unit cells structures required to compute the leads self-energies (see below). All the structures were relaxed through a standard conjugate gradient algorithm until all the forces were smaller than 0.04 eV/Å.

In the transport calculations we partition the system into three regions: a left lead, a right lead and a central scattering region that contains the twisted section of the ribbon. We solve the electronic transport problem in the central region with the transiesta method [55] within the nonequilibrium Green’s function formalism. The open boundary conditions imposed by the electrodes are accounted for through the left (right) self-energy $\Sigma_{L,R}(E)$. The zero-bias transmission $T(E)$ is calculated as

[TABLE]

where $\mathbf{G}_{C}^{r,a}(E)=[E\mathbf{S}_{C}-\mathbf{H}_{C}-\mathbf{\Sigma}_{L}^{r,a}(E)-\mathbf{\Sigma}_{R}^{r,a}(E)]$ is the retarded (advanced) Green’s function of the scattering region, $\mathbf{\Gamma}_{L,R}(E)=i(\mathbf{\Sigma}_{L,R}^{r}-\mathbf{\Sigma}_{L,R}^{a})$, $\mathbf{H}_{C}$ is the Hamiltonian matrix and $\mathbf{S}_{C}$ is the overlap matrix. We use three unit cells of the armchair GNR for both electrodes (see Fig. 1 (a)), while the extent of the scattering region varies and depends on the torsional deformation considered in each case. The electrical conductance $G^{e}(\mu)$ is then calculated through the Landauer formula as

[TABLE]

where $f_{0}(E,\mu)$ is the Fermi-Dirac distribution function at a chemical potential $\mu$.

2.2 Phonon transport

We simulate phonon transport with nonequilibrium Green’s function methodologies on the same lead/scattering region/lead structures used for the electron transport calculations. For phonons, however, the use of DFT would be impractical due to the large number of calculations needed to compute the force-constants matrix. Instead, we relax the structures and calculate the forces using the classical bond-order potential due to Brenner [56] as implemented in the GULP code. [57]

The phonon transmission function, $\mathcal{T}^{ph}(\omega)$, is obtained from the phonon equivalent of Eq. 1, where now the retarded (advanced) Green’s function of the scattering region reads $\mathbf{G}_{C}^{r,a}(E)=[\omega^{2}\mathbf{I}-\mathbf{F}_{C}-\mathbf{\Sigma}_{L}^{r,a}(\omega)-\mathbf{\Sigma}_{R}^{r,a}(\omega)]$, with $\mathbf{F}_{C}$ being the force-constant matrix, $\mathbf{I}$ the identity matrix and $\mathbf{\Sigma}_{L,R}(\omega)$ the self-energy of the left (right) contact. In analogy with Eq. 2, the phononic thermal conductance is calculated within Landauer theory as

[TABLE]

where $n_{0}$ is the equilibrium Bose-Einstein distribution function.

Besides, we employ a recently proposed [58] extension of the nonequilibrium Green’s function method to disentangle the contribution of each individual phonon mode to the transmission function and the thermal conductance. To this end we calculate the Bloch matrices of the leads as explained in Ref. [58] whose eigenvalues and eigenvectors provide access, on the one hand, to the dispersion relations and eigendisplacements of the phonons at both leads and, on the other hand, to a single-mode transmission matrix, $\mathbf{t}(\omega)$. The square modulus of an element of this transmission matrix, $|t_{m,n}(\omega)|^{2}$, represents the probability of transmission from the $m$th phonon mode in the right lead to the $n$th mode in the left lead at a given frequency and has a value between 0 and 1. The sum of all single mode transmissions equals the result given by Eq. 1.

Our transport simulations assume a ballistic regime, i.e., they are robust as long as electron-electron, phonon-phonon and electron-phonon interactions can be neglected. Despite such limitations, this a reasonably clean approach where the only source of scattering is the structural deformation of the ribbon.

2.3 Thermoelectric properties

The thermoelectric figure of merit can be defined as,

[TABLE]

where $S$ is the Seebeck coefficient, $G^{e}$ is the electronic conductance, $T$ is the temperature, and $G^{ph}$ and $\kappa^{e}$ are the thermal conductance due to phonons and electrons, respectively. We calculate $G^{e}$ and $G^{ph}$ as indicated in Eqs. 2 and 3, whereas $S$ and $\kappa^{e}$ are calculated from the electron transmission function (Eq. 1) as follows [59, 60],

[TABLE]

[TABLE]

where the functions $L_{m}(\mu)$ are defined as,

[TABLE]

3 Results

3.1 Electron transport

We start our study by considering an armchair GNR with a full twist of $180^{\circ}$ that develops over increasing lengths, $L_{twist}$. Our goal is two-fold: (i) estimating the minimum length required so that a $180^{\circ}$ twist does not yield a significant decrease of the conductance (in the limit of an infinitely long scattering region the behaviour of the twisted system must approach that of the flat GNR); (ii) characterizing the degradation of the conductance at those short lengths where the effect of the distortion is not negligible.

Some examples of the systems studied are sketched in Figure 1. We have analysed 12 different lengths of the twisted region, ranging from $L_{twist}=9$ to 20 unit cells of the armchair GNR, i.e. from 38.5 to 85 Å. Notice that, prior to the calculation of the self-energies of the leads that allow treating the system as semi-infinite, the electronic structure within periodic boundary conditions must be calculated [55]. However, if the GNR is symmetric with respect to the twist axis, a full twist of $180^{\circ}$ naturally satisfies the requirement of periodicity. This is not the case with the arbitrary twist angles that will be studied in the second part of this section.

The zero-bias conductance of a subset of the investigated systems is shown in Fig. 2 together with that of the undistorted ribbon. We observe that the shorter is the length within which the full twist is forced to occur, the lower the conductance, with the most evident changes occurring for energy values next to the band edges and to the increase from one to two transport channels ($\sim 1$ eV above (below) the bottom (top) of the conduction (valence) band). On the other hand, large values of the central length determine a generally high conductance, eventually approaching the value of the undistorted ribbon in the asymptotic case of infinite length. As matter of fact, we obtain very similar conductances for lengths of the twisted regions larger than 14 unit cells (see, e.g., the cases with $L_{twist}=$ 15, 17, and 19 in Fig. 2).

Additional insight of the scattering mechanism is given by the analysis of the transport eigenchannels. We have plotted the eigenchannel that carries most of the current at two selected energies for a short and a long GNR (Fig. 3). At $E=0.2$ eV, the short GNR suffers a considerable scattering and very little of the incoming eigenstate from the left-hand side survives and is found on the right-hand side, after the $\pi$ twist. For the longer GNR, on the other hand, the eigenstate is transmitted, although after crossing the twist it changes its symmetry. At $E=0.47$ eV the transmission is very high for all the GNR investigated, as we already know from the analysis of the transmission. Accordingly, the eigenchannel transmits well and similarly for both considered cases, maintaining the same symmetry and a similar magnitude.

In view of the results discussed above, we now give a closer look at what twisting a GNR means from the structural viewpoint. When a nanoribbon is twisted, it undergoes two main modifications with respect to the planar configuration. On one hand, the structure is no longer flat; on the other side, the torsional distortion alters the bond-length distribution. Specifically, in an armchair GNR, twisting affects mainly the components of the interatomic distances parallel to the transport direction, i.e. the ribbon axis. The extent to which the bond-lengths are modified is related to the length of the region over which the twisting takes place: the shorter the nanoribbon, the stronger the impact on the distances among the atoms. The dependence of this structural distortion on $L_{twist}$ is shown in fig. 4, where the change in the bond length with respect to the flat ribbon is plotted as a function of the longitudinal position. The value is averaged over the number of atomic bonds in a given spatial interval. Each histogram corresponds to a specific length of twisted nanoribbon: a short, a medium and a long system ($L_{twist}=11$, 16 and 20, respectively). The averaged deviation in the interatomic distance is almost uniform along the length of the twisted region (while it expectedly goes to zero at the extremes of the central scattering region where the geometry is kept frozen during the relaxation to provide the correct coupling with the leads in the transport calculations). More importantly, however, Fig. 4 shows how the variation in interatomic distances tends to decrease for increasing lengths. A more concise representation of these results is given in the inset, where by plotting the standard deviation of the bond length, we quantify the amount of local distortion as a function of $L_{twist}$.

We now study the dependence of the conductance of an armchair GNR in presence of a torsional deformation of an arbitrary angle. As mentioned above, before coupling the scattering region to the leads a calculation of its electronic structure in periodic boundary conditions is needed. Notice that such a requirement of periodicity of the twisted system is only guaranteed in the case of a torsion of $180^{\circ}$. To avoid this limitation we study GNRs where a twist of $\phi$ that develops in the first half of the scattering region is followed by an opposite twist of $-\phi$ in the remaining half. Some atomistic structures of the studied systems are shown in Fig. 5. Notice that, at variance with the case of a full twist, here we need to freeze a small region of the GNR, right where the twist direction is reversed, to prevent the system to relax back to the flat ground state. We consider a fixed length of 16 unit cells of the scattering region and vary the twist angle.

Before performing a systematic study of the dependence on the angle, we have analysed which is the effect of a reversal of the twist angle, i.e. a torsion of $\phi$ followed by a torsion of $-\phi$. In particular, we have considered the case of $\phi=90^{\circ}$: a twist of $90^{\circ}$ in the first half of the GNR, followed by a twist of $-90^{\circ}$. We compare the transmission function of this case with a single, continuous twist of 180∘ that develops all over the length of the scattering region. The results are shown in Fig. 6 (a). As can be seen there, the two curves are almost indistinguishable. This means that reversing the twist angle has virtually no effect and a series of twists of $+\phi/-\phi$ is roughly equivalent to a continuous twist of $2\phi$. Therefore, we will present the results that follow as a function of $\theta=2\phi$, allowing direct comparison with the results discussed in the first part of this section.

Fig. 6 (b) shows the conductance of the twisted GNR for different twisting angles $\theta$ in the range $[0,360^{\circ}]$. As can be seen there, an increase of $\theta$ determines an increase of scattering in the central region with a following reduction in conductance. The curves at the highest twist angles ($\theta\geq 240^{\circ}$) exhibit a strong chemical potential dependence, with a well recognizable peak structure that indicates that only at specific energies the transmission is finite, whereas it is almost fully suppressed elsewhere.

Therefore, it seems that for sufficiently large twist angles there is a transition from a band-like one-dimensional conductance, to a transport mechanism typical of molecular systems or constrictions, where energy is transmitted only in correspondence of resonant levels. To better look into this transition, once again we study the transmission eigenchannels and their spatial distribution, which we plot in Fig. 7 at two selected energies. At $E=0.27$ eV, where the transmission is high in all cases (see Fig. 2) no difference can be appreciated in the cases considered ($\theta=160^{\circ}$, 280∘, and 360∘). The eigenchannels are mostly delocalized over the twisted region and guarantee a good coupling with the electrode incoming states, thus resulting in a high transmission. On the other hand, moving to $E=0.42$ eV, where the conductance of GNR with high twist angles has a minimum, the eigenchannel is extended for small angles but becomes more and more localized as $\theta$ grows. As it can be seen, there is an increasingly high probability of finding the electrons in the region between two knots, where the change in the rotation direction takes place. This localized electronic states still allow electron transport, but only at resonant energies with these knot levels.

3.2 Phonon transport and thermoelectric figure of merit

The effect of twisting a graphene nanoribbon on its phonon-mediated thermal transport properties has been previously studied with molecular dynamics [61, 62] and NEGF simulations. [63] These studies have shown that, as for the case of electrical transport, the thermal conductance gets gradually reduced as a stronger torsional deformation is applied on the nanoribbon. In addition, the following trends have been observed: i) The conductance of zigzag terminated ribbons is more sensitive to twisting than that of armchair ones, [62] ii) the conductance is more effective reduced by twisting at low temperature regimes, when conductance degradation due to anharmonic phonon-phonon scattering is less relevant, [61] iii) twisting a nanoribbon has a sharper effect for narrower and shorter ribbons. [62, 63] In spite of these previous studies on the topic, it remains as an open question to determine the degree in which the torsional deformation affects the transport of specific phonon polarizations.

Here we perform mode-dependent transport simulations and calculate the contribution to the thermal conductance of twisted GNRs from longitudinal (L), transverse (T) and out-of-plane (Z) phonons. To classify a phonon-mode as L, T or Z we adopt as a practical criterion that its eigenvector presents at least 70% polarization along one of the three directions. In Fig. 8 we show the results of the polarization-resolved conductance analysis as the twisting angle is increased for a fixed nanoribbon length of 16 unit cells (upper panels) and as the length of the twisted region is modified for a fixed twisting angle of $180^{\circ}$ (lower panels). The corresponding calculation on a flat nanoribbon yields that Z and L phonons contribute much more than T phonons to the thermal conductance (not shown). Upon twisting the ribbon we observe that the conductance due to Z phonons experiences the largest degradation, of the order of 70% at 360*∘*, through a process that shows large independence on the temperature. This contrasts with the results of L and T phonons for which twisting the ribbon reduces more effectively their contribution to the conductance at higher temperatures. On the other hand, upon shortening the twisting region we observe that the conductance due to L phonons is the most affected with respect to the flat nanoribbon case.

Finally we exploit the obtained results to assess the performance of twisted GNR systems for thermoelectric applications. In this field, a search for strategies to enhance charge transport and/or deteriorate heat conductivity is currently ongoing for a variety of materials.  [64, 65, 66, 67, 68] A measure for the performance of thermoelectric systems is given by the figure of merit ZT (Eq. 4). In this respect, it is clear that the relative change of electron and phonon conductance caused by twisting in our graphene nanoribbons can be responsible for a large modification of thermoelectric figure of merit with respect to the undistorted system. Indeed, depending on the relative variation of the different ingredients of ZT, an improvement or deterioration in the performance will follow.

The evolution of ZT with the chemical potential at T=50 K is shown in Fig. 9 for different twisting lengths and angles. Two general trends are observed. First, the chemical potential can be tuned in order to increase the value of ZT, e.g. by applying a gate voltage to the twisted channel. Second, the ZT is shown to have peaks in correspondence of conduction and valence band edges. [69] As can be seen in Fig. 9 (a), the ZT peaks decrease when the twisting length is reduced and their value is always smaller than that of the undistorted ribbon. The figure of merit is expected to tend asymptotically to this maximum value for sufficiently large twisting lengths. Regarding the effect of the twist angle (Fig. 9 (b)), we observe that ZT peaks are also found at energies close to the band edges for sufficiently small angles. However, when the twisting is made stronger ($>280^{\circ}$), a new collection of peaks appears at new values of chemical potential corresponding to maxima of electronic conductance (see Fig. 2). This behaviour consistently reflects the discussed change in the nature of electron transport in these systems, evolving from 1D semiconducting channels to molecule-like systems. Notwithstanding, it is important to note that, also for this differently twisted ribbons, the figure of merit is smaller than in the untwisted case. We have checked that this general observation is ruled by the reduction of the Seebeck coefficient as the GNR is deformed.

4 Conclusions

In this work we simulated ballistic electron and phonon transport in graphene nanoribbons with torsional deformations through a nonequilibrium Green’s functions procedure and determined the following. Broadly speaking, both electrical and thermal conductances are reduced by twisting a nanoribbon. The shorter is the length within which a given twist is applied, the larger is the structural deformation, which is reflected in the distribution of atomic bonds lengths, and the lower is the electrical conductance. In turn, deformations applied over a long portion of nanoribbon have a negligible effect on the transport properties. On the other hand, when two consecutive twists are applied in opposite angular directions a transition from a 1D to a 0D (molecular) electron transport regime is observed. The molecular-like transport is a consequence of the formation of localized states in the region of the nanoribbon between the two knots. Concerning phonon transport, we have observed that the propagation of out-of-plane and longitudinal vibrational modes is the most impeded by increasing the twist angle and by shortening the twist length, respectively. Finally, the thermoelectric efficiency of GNRs is generally attenuated by torsional deformations because these entail a reduction in the Seebeck coefficient, thereby reducing the figure-of-merit ZT.

Acknowledgements

We thank Troels Markussen for useful discussions. M.R. and R.R. acknowledge financial support by the Ministerio de Economía y Competitividad (MINECO) under grant FEDER-MAT2013-40581-P and by the Severo Ochoa Centres of Excellence Program under Grant SEV-2015-0496 and by the Generalitat de Catalunya under grants no. 2014 SGR 301 and through the Beatriu de Pinós fellowship program (2014 BP_B 00101). We thank the Centro de Supercomputación de Galicia (CESGA) for the use of their computational resources.