Modelling of energy transfer induced by longitudinal shaking of one atom row in carbon nanoribbons

1Institute of Molecule and Crystal Physics URC RAS, Prt. Oktyabrya 151, 450075, Ufa, Russia 2Institute for Metals Superplasticity Problems, RAS, Khalturina 39, 450001, Ufa, Russia 3Department of Physics, Ferdowsi University of Mashhad, Vakilabad Highway, 9177948974, Mashhad, Iran 4Bashkir State Agrarian University, 50 let Oktyabrya 34, 450001, Ufa, Russia 5Peter the Great St. Petersburg Polytechnic University, 29 Polytechnicheskaya St., 195251, St. Petersburg, Russia


Introduction
The concept of discrete breathers (DB) have become very useful for describing localized excitations in defect free nonlinear lattices. Starting from the work of Sivers and Takeno [1] nonlinear dynamical objects of that type drew significant attention of the scientific community. Discrete breathers (also referred to as intrinsic localized modes or ILM) were directly observed in many macroscopic and mesoscopic system of different physical nature [2][3][4][5][6]. On the other hand, the theoretical investigations of discrete breathers have been also developed to demonstrate that in general every discrete nonlinear system with translational symmetry can maintain breathers [2,7,8]. The frequency of DB lies outside of the system phonon spectrum, which means that once excited it will not emit energy in the form of plane waves. Nevertheless the very fine tuning of the form of initial excitation is needed to introduce a perfect breather (which will not radiate energy) to the system [9,10]. In the physical experiments it is only possible to have breather-like localized object which can be referred to as quasi-breather [10,11] and has a finite lifetime because of the energy loses during the oscillations.
The fact that discrete breathers have very weak interactions with phonon modes makes them very promising objects to be used in crystals [12]. Despite of there are several works on experimental study of DBs in crystals [13][14][15][16][17], their direct observation in such systems is complicated by the presence of thermally activated delocalized oscillations [18]. That enlightens computer modelling of breathers in crystals as one of the most preferable approach of their study. By means of molecular dynamic (MD) the properties of discrete breathers have already been established in several crystal lattices [19][20][21][22][23][24][25] Such calculations are often held in the classical framework based on the phenomenological interaction potentials. The choice of particular potential is a very sophisticated problem especially for the modelling of highly nonlinear dynamical objects such as DB. To give an example of differences in the breathers properties due to the use of different phenomenological potentials we refer to the work [26] where authors investigate the intrinsic localized modes in Si crystal. On the other hand in the work [8] it was shown that properties of discrete breathers in graphane (fully hydroginated graphene) obtained from the ab initio calculation coincide with those obtained in the classical MD [27] only in the small range of low DB amplitudes. The latter results were received with the use of the well known AIREBO potential [28].
In graphene the existence, stability and properties of discrete breathers have been studied both by means of classical MD [29][30][31] and ab initio calculations [32]. The latter approach overcomes the problem of phenomenological potential choice and can be used to verify the results of the former approach. It is shown in [32] that the use of the standard set of potentials for carbon and hydrocarbon systems [33] leads to the properties of DBs which are in a good agreement with those obtained from first principal calculations. That fact gives an opportunity of studying highly nonlinear excitations in graphene by means of classical molecular dynamics on relatively long time intervals that could not be achieved in ab initio calculations. Moreover, thermoconductivity of graphene nanoribbons is still not investigated well but already is considered as very high, which means that a thorough study of the energy flow in graphene is needed.
The present paper is devoted to the investigation of energy transfer mechanism in graphene nanoribbons by means of classical molecular dynamics simulations. The calculations are held using the standard set of interatomic interaction potentials in carbon and hydrocarbon structures [33] (also sometimes referred to as Savin potential). A distinctive feature of the potentials is that they reproduce the dispersion curve of graphene better than the Brenner potentials [34].

The model
A graphene sheet forms a two-dimensional hexagonal crystal lattice with a primitive cell containing two atoms. To simulate the nanoribbons, a rectangular supercell with dimensions M×N was built, where M and N are number of unit cells along x ("zigzag" graphene direction) and y ("armchair" graphene direction) axes, respectively. The periodic boundary conditions were applied to the calculation block (see Fig. 1a). The standard set of interatomic potentials gives the valent bond length in the unstrained graphene ρ 0 = 1.418 Å. Without any strain applied the phonon spectrum of graphene doesn't contain a gap. As long as previous studies [29,32] showed that discrete breathers in graphene are dynamical objects with soft type nonlinearity (their frequency is decreasing with increase of amplitude), we introduce a strain with the components ε xx = 0.35, ε yy = −0.1, and ε xy = 0.0. The equilibrium positions of atoms in uniformly strained graphene are found by minimizing the potential energy of the crystal. For the chosen strain components the equilibrium flat configuration of graphene is stable [35]. The phonon density of states of our system under the strain is depicted in Fig. 1b. It shows a wide gap in xy phonon states which makes it possible for DBs to exist with frequencies inside the gap (see [29] for the thorough discussion).
To start an energy transfer to the system we introduce a shaking of the chain of the carbon atoms in the middle of the ribbon. Since our system is longer on the y direction then in x direction, the atoms which are being shaking form a zigzag shaped chain. The number of atoms in that chain equals 2M. We choose now the most simple type of the shaking, namely the "displacement driving" one. For that case the exact displacements of the shaking atoms are in the form ∆y = Asin(2πνt), where ∆y represents the displacement of atoms from their equilibrium positions in y direction, A is the amplitude in angstroms and ν is the frequency in THz. The main quantity which describes the energy transfer is power P = ∆E/∆t, where ∆E stands for the additional energy which came to the system through the time ∆t.

Results
For the displacement driven shaking we study different frequencies ν and amplitudes A of the driving. For each pair of parameters the time of calculation was estimated to let the excitation go just to the last particle on the each end of the graphene ribbon. Further calculations would be affected by the border effects. In the Fig. 2  The dependence of power on frequency (Fig. 2) should be compared with the phonon density of states (Fig. 1b). It is obvious that the main regions, where the energy transfer to the system is possible, coincide with the regions of xy phonon waves existence (frequencies from 0 to 22.5 THz and from 35.4 to 40.3 THz). Nevertheless there is the region where frequency of shaking doesn't match with any phonon wave frequency (in the gap of phonon spectrum), but still energy transfer is possible with not negligible power values. The inside picture on Fig. 2 shows the magnification of the frequency region for which shaking frequency lies in the band gap of graphene phonon spectrum. For that case energy goes to the lattice in a very special manner. On the Fig. 3 the evolution in time of two particles closest to the shaking chain is shown for two different frequencies. The Fig. 3a corresponds to the driving frequency ν = 26.7 THz. It can be seen that the maximum amplitude of the atom closest to the shaking chain is 10 times higher than the driving amplitude. Nevertheless, the breather continuously loses and gains energy through the interaction with the shaking row. The result on the Fig. 3a does not coincide with the sharp peak of power transfer, located at ν c = 26.9 THz (see inside picture on Fig. 2). The energy transfer for ν = 26.7 THz is very small.
The analysis of the oscillation frequencies of atoms close to the shaking chain shows that for the particular amplitude A of the shaking the critical frequency ν c corresponds to the frequency of discrete breather which appears in the system. For the case of shaking parameters A = 0.01 Å and ν = 26.7 THz the frequency of DB floats close to the mean value of ν db = 27.3 THz with standard deviation σ νdb = 0.1 THz (see oscillations on Fig. 3a). For the case of shaking parameters A = 0.01 Å and ν = 27.7 THz (above the ν c = 27.2 THz for that amplitude) the frequency of DB floats close to the mean value of ν db = 26.7 THz with standard deviation σ νdb = 0.6 THz (see oscillations on Fig. 3b). As it can be seen, if the frequency of shaking lays below the frequency of DB which appears in the system, the energy transfer to the nanoribbon is very low, so the breather prevents the energy flow. In the case of driving frequency above the ν db , breather promotes the energy spreading.

Conclusion
The simple case of displacement driven energy transfer to the carbon nanoribbon is investigated. It is shown that energy can spread over the system for shaking frequencies not only within phonon spectrum of graphene. Inside the gap of phonon spectrum there is a resonant frequency ν c for which the power of energy transfer reaches highest values. The frequency ν c coincides with the frequency of the DB emerging close to the shaking chain. For the given amplitude of shaking if the driving frequency is less then the DB frequency, than energy transfer to the system is very low. In the case of shaking frequencies higher than ν c the discrete breather promotes energy flow to the ribbon due to the low frequency modulations of its' amplitude.