Asymmetric energy transport in defected boron nitride nanoribbons : Implications for thermal rectification

Related Articles Intrinsic nanoscale phase separation in miscible mixtures and supramolecular assemblies AIP Advances 1, 042173 (2011) Ultra-high four wave mixing efficiency in slot waveguides with silicon nanocrystals Appl. Phys. Lett. 99, 191105 (2011) Limits of elemental contrast by low energy electron point source holography J. Appl. Phys. 110, 094305 (2011) Nitrogen flux induced GaN nanostructure nucleation at misfit dislocations on Al2O3(0001) Appl. Phys. Lett. 99, 143105 (2011) Core size dependence of the confinement energies, barrier heights, and hole lifetimes in Ge-core/Si-shell nanocrystals J. Appl. Phys. 110, 074306 (2011)


I. INTRODUCTION
5][6][7][8] While the electronic properties of nano-BN structures have been well characterized, it is only very recently that the phonon/thermal properties have been systematically investigated, primarily using computations as the probing tool.0][11][12][13] Of particular interest is the work by Yang et al., 13 who used non-equilibrium Green's function methods to examine the phonon transport properties of armchair-edged BNNR (a-BNNR) and zig-zag BBNR (z-BNNR) that contained triangular vacancies, first observed by Jin et al. 14 in electron-irradiated SLBN.A key result was the fact that the thermal conductance of BNNRs decreased linearly with increasing size of the triangular vacancy; equally importantly, as pointed out by Yang et al., this opens up new avenues for manipulating the thermal properties of defected BNNR.In particular, the geometric asymmetry of the triangular vacancies can be exploited to preferentially scatter phonons in BNNR, which can, in turn, lead to spatially-dependent thermal properties as dictated by the orientation of the triangular vacancies.A related consequence is the possibility of realizing thermal rectifiers as discussed by Go et al., 15 where it was shown that the thermal conductivity of a material has to be an inseparable function of both space and temperature to exhibit thermal rectification.In this context, developing a fundamental understanding of the thermal transport properties as a function of the orientation of the triangular vacancy forms the basis of this work.Towards this end, we employ molecular dynamics (MD) simulations to examine the interplay between defect-orientation and thermal-energy transport properties of BNNR.

II. MODELS AND METHODS
The accuracy of MD simulations is governed by the choice of the underlying interatomic potentials and Brenner-Tersoff-style potentials have been successfully used to model BN systems, some which have been parameterized to accurately yield the phonon dispersion of SLBN. 9,11 n this study, the Brenner-Tersoff parameters as developed by Albe and Moller 16 will be used as they have been successful in modeling the different hybridization-states of BN-an important requirement while modeling defected BNNR.Formally, the Brenner-Tersoff potential is given by where, E i is the energy of atom i, V ij is the interaction energy between atoms i and j, and n i is the number of 'interacting' neighbors of atom i.The parameters are given in Ref (16).The primary materials-parameters that define thermal transport are thermal diffusivity and thermal conductivity.Strategies to evaluate these parameters include non-equilibrium (NEMD) and equilibrium (EMD) methods.In the NEMD framework, the thermal conductivity is obtained directly by solving Fourier's law under steady-state conditions, where a temperature gradient is maintained across the modeled material by fixing the temperature of the two ends of the material at different temperatures.Thermal diffusivity is evaluated under transient conditions, by solving for the second-order heat equation.The EMD method is based on the Green-Kubo formulations, where NVE (i.e. the microcanonical ensemble) conditions are imposed on the simulated system; based on the equilibrium fluctuations in the heat current (S), the thermal conductivity (κ) is estimated from the time-dependent autocorrelation of the heat current (HCAF) as given by Eqn.(2), where V and T are the volume and temperature of the system respectively. 17The thermal diffusivity (D T ) can then be calculated (see Eqn. ( 3)), where C p and ρ are the specific heat and mass density respectively.
EMD and NEMD methods have been routinely used to model thermal transport in materials, but care has to be taken in their implementation; NEMD methods imposes extraordinarily large temperature gradients across the material that may not be realized experimentally; further, as discussed by Jiang et al. 18 as also observed by the authors, 19 the thermostated ends induce spurious vibrational modes characteristic of the size and location of the respective thermostats, which modify the injected heat flux, leading to the possible erroneous estimation of the thermal conductivity.While EMD can yield the thermal conductivity of materials based on the Green-Kubo methods, it is not straight-forward to determine direction-dependent conductivity, which is necessary to characterize the effect of triangular-vacancy orientation on thermal transport.Also, while modeling dimensionally-restricted systems such as BNNR, where the phonon mean free path can be comparable to the system-size (i.e.ballistic phonon transport), the definition of temperature becomes debatable. 20In addition, MD simulations use the classical equipartition definition of temperature, which is strictly true only when above the Debye temperature, where the Bose-Einstein distribution reduces to the classical Boltzmann definition. 20Clearly, in order to obtain quantitatively accurate estimates of thermal conductivity and diffusivity of BNNR, the above shortcomings have to be unambiguously addressed; alternately, if appropriate interpretation of MD results is carried out, the interplay between defect orientation and thermal transport in BNNR can be suitably examined without explicitly quantifying conductivity and diffusivity of pristine and defected BNNR.Towards this end, EMD simulations will be employed to study the thermal conductivity characteristics.Specifically, by calculating spatial variations in the heat current autocorrelation function (HCAF), a relative measure of the spatial variation in thermal conductivity and its interplay with defect-orientation can be obtained, leading to verifying the possibility of defected BNNR exhibiting thermal rectification characteristics.In order to study thermal diffusivity, a variant of NEMD simulations will be undertaken.Here, one end of the BNNR is suddenly quenched and held at a fixed temperature; the time taken for the temperature of the rest of the material to equal that of the thermostated end is taken as a measure of the thermal diffusivity.By using this method, we avoid estimating heat fluxes, and thereby circumvent associated problems discussed above.Thus, using the above described methods, indirect measures of the steady-state (i.e.thermal conductivity) as well as the transient thermal properties will be examined, enabling the characterization of the effect of extended triangular vacancies on thermal transport of pristine and defected BNNR.It should be noted that in this work we use the equipartition definition of temperature and do not use quantum corrections, given that we only rely on qualitative and not quantitative comparative analysis to deduce the thermal transport behavior.

III. SIMULATION PROCEDURE
The length of the simulated a-BNNR was 17.5 nm, while periodic boundary conditions were applied along its 7.1 nm width.Fixed boundary conditions were imposed on the edge atoms (i.e. the thinnest strip consisting of boron and nitrogen atoms at each end).For the defected system, the defect was represented by an equilateral triangle with nitrogen-termination to ensure consistency with experimental observations.The defect orientation is shown in Fig. 1, and its dimensions were chosen to be approximately half the BNNR width (corresponding to a filling fraction of 5%).
In order to carry out EMD simulations, both pristine as well as defected BNNR were initially equilibrated under NPT conditions at 300 K and 900 K respectively, after which NVE simulations were conducted for 7 ns with a time step of 1 fs.The systems were spatially divided into 15 bins along their length to enable the calculations of spatial variations in HCAF, which were obtained over the last 5 ns of each NVE run.Particular attention was paid to the HCAF component along the length of the BNNR, which was primarily used in our data analysis.To ensure better statistics, five different equilibrated starting configurations were used for each case.In the NEMD simulations, the thinnest possible strip of atoms (consisting of equal number of boron and nitrogen), adjacent to the boundary atoms at the opposite ends of the BNNR, were identified to be the thermostated regions which were governed by a Nose-Hoover thermostat.The boundary atoms were not included to avoid edge effects as noted by Jiang et al.. 18 For the 300 K and 900 K systems, the thermostat temperature equaled 150 K and 450 K respectively, and the time for the rest of the unconstrained system to attain the temperature of the thermostated region was calculated when the thermostat was placed at the (i) left and the (ii) right edge respectively.In particular, consider the first two HCAF valleys/peaks in the bins between the triangle-face AB and the near edge (i.e.bins 1-6).The bins in proximity to the BNNR edge (bins 1-3) are characterized by valleys followed by peaks in HCAF, while the HCAFs in bins closer to the triangleface AB (bins 4-6) are first described by peaks and subsequently by valleys.This is explained by the fact that the triangle-face AB is not a fixed-boundary, and does not lead to a phase-change during reflection.At approximately 0.7 ps, we see the uniform appearance of valleys for all the bins (1-6), which is correlated to the simultaneous 'round-trip' arrival of two phonons.Note that the same phenomenon is also observed in the pristine BNNR system at 1.5 ps-due to the longer path (almost twice) traversed by the respective phonons.For the region in the defected BNNR between triangle vertex-C and the farther edge, the spatial HCAF profile diverges from that of the other region; all the bins corresponding to this region (9-15) are characterized by an initial valley (reflection from the farther fixed end), but subsequent features are not well pronounced, a direct consequence of phonon scattering from the sloped edges of the triangle defect, which can be distinguished from the reflection that occurs at the normal AB face.Thus, phonon propagation characteristics in the two regions separated by the geometrically-asymmetric triangular defect are indeed different.
Figs. 4(a) and 4(b) illustrate the HCAF of pristine and defected BNNR at 900 K.A comparison with Figs. 2 and 3 indicates the role of temperature on the HCAF profile.Clearly, the anharmonic effects become more distinct at the higher temperature, as seen by the absence of higher order HCAF echoes in the respective systems.Thus, by comparing and contrasting the HCAF characteristics of pristine and defected BNNR, one can conclude that geometric asymmetry of the defect leads to distinct spatial-and temperature-dependent thermal-phonon propagation characteristics for the defected BNNR system, indicating the possibility of observing thermal rectification in such systems.
In order to study the transient response of the two systems, the quenching procedure as described earlier was adopted.Fig. 5(a) and Fig. 5(b) illustrates the rate of temperature change for the pristine BNNR and defected BNNR (both initially at 300 K and quenched to 150 K), when in 'forward' bias (i.e. the thermostated BNNR edge faces the triangle-face AB) and 'reverse' bias (i.e. the thermostated BNNR edge faces the triangle vertex-C).While the pristine BNNR responds identically under both forward and reverse bias, the temperature-time curve do not overlap for the defected BNNR, implying that the thermal diffusivity is position-dependent.Numerical solution of the transient heat equation shows that the reverse-bias apparent thermal diffusivity is higher by a factor of 1.13.A similar result was also observed when the 900 K systems were quenched to 450 K, with the ratio of the reversebias to forward-bias thermal diffusivity for defected BNNR equaling 1.07.These results when viewed in conjunction with the HCAF observations clearly indicate that the asymmetric triangular defect plays an important part in the ability of the defected BNNR to respond to external thermal stimuli.Specifically, based on the orientation, specific triangular-faces can impede phonon-energy propagation, thereby allowing defected BNNR systems to exhibit spatially asymmetric thermal transport properties.

V. CONCLUSIONS
Using robust atomistic models, we have examined the effect of geometrically-asymmetric triangular defects on thermal transport in BNNR.By suitably interpreting the time-evolution of spatially decomposed HCAF in terms of phonon-propagation, we have unequivocally demonstrated that the triangular defects can significantly affect the thermal transport characteristics.This was further confirmed by appropriate analysis of direction dependent thermal diffusivity estimations in BNNR.This study signifies an important first step in enabling the use of BN-based materials for important thermal management applications such as thermal interface materials and thermal rectifiers.Nevertheless, comprehensive and detailed studies have to be undertaken to completely characterize the interplay between the BNNR dimensions, triangular-defect size and orientation and the thermal properties.In particular it is important to identify the critical BNNR dimensions where intrinsic Normal and Umklapp processes become dominant (as compared to defect-phonon scattering) leading to isotropic, diffusive thermal transport.Such studies can enable the experimental realization of BN systems with tunable thermal properties.

FIG. 2 .
FIG. 2. (left)-Time evolution of the spatially-resolved HCAF for pristine BNNR at 300 K. Black overlays are an aid to the eye.(right)-Representation of the spatial decomposition of the simulated system into 15 bins.

FIG. 3 .
FIG. 3. (left)-Time evolution of the spatially-resolved HCAF for defected BNNR at 300 K. (right)-Representation of the spatial decomposition of the simulated system into 15 bins.