A Semi-Grand Canonical Kinetic Monte Carlo study of Single-Walled Carbon Nanotubes growth

Single-walled carbon nanotubes exist in a variety of different geometries, so-called chiralities, that define their electronic properties. Chiral selectivity has been reported in catalytic chemical vapor deposition synthesis experiments, but the underlying mechanisms remain poorly understood. In this contribution, we establish a simple model for the prediction of the growth rates of carbon nanotubes of different chiralities, as a function of energies characterizing the carbon nanotube-catalyst interface and of parameters of the synthesis. The model is sampled efficiently using kinetic Monte Carlo simulations in the semi-grand canonical ensemble, uncovering the interplay of the external experimental conditions and the configuration and energetics of the interface with the catalyst. In particular, the distribution of chiral angle dependent growth rates follows non-monotonic trends as a function of interface energies. We analyze this behavior and use it to identify conditions that lead to high selectivity for a variety of chiral angles.


I. INTRODUCTION
Single-Walled Carbon Nanotubes (SWNTs) are crystalline structures that can be seen as rolled up graphene sheets with a remarkably large aspect ratio, reaching up to tens of centimeters in length and a few nanometers in diameter.Since their first observation and characterization 1 tremendous research efforts went into the development of efficient synthesis techniques, the most prominent being the catalytic chemical vapor deposition (CVD), in which a carbon-rich gaseous feedstock (e.g., CO, CH 4 , ethanol) is being decomposed at high temperature, on the surface of a catalyst, often a transition metal, to nucleate and grow a tube 2 .A number of groups [3][4][5][6] have reported the selective growth of tubes with a particular helicity (or chirality) characterized by the so-called Hamada indices (n, m).However, we are far from understanding all aspects of this control of chirality which, if routinely achieved, would promote the development of technologies based on the exceptional electronic properties of SWNTs 7 .A reason for the difficulties with selective growth is that SWNTs may appear in a (theoretically) infinite number of polymorphs with very small energy differences between them, pushing the crystal growth science to its limits.
While SWNTs are indeed crystals, they display characteristic properties making them different from 3-dimensional materials.The first obvious difference is that the active area for carbon incorporation in the tube during CVD growth is made of only a few tens of atoms at the interface between the tube and its seeding catalyst.Because of the high temperature of the CVD synthesis and of its nanometric size, this interface is experimentally difficult to characterize and even more to control.A second aspect is that some local topological defects, such as a pair of non-hexagonal carbon rings, have the potential to change the crystal structure of the tube 8 .Point defects are present and thermodynamically stable in 3-dimensional materials without altering their crystallographic structure.However, this is not always the case in the present instance, and the fact that perfect SWNTs with a defined structure can be grown up to centimeter length 9,10 , implies a negligible concentration of topologically active defects in this size range.
Classically, crystals are formed by nucleation and growth.The former is dominated by the thermodynamic stability of the nucleus, while kinetics is central to the latter.We recently developed a thermodynamic model of the catalyst-tube interface where the fluctuations of the interface structure and the associated entropy were found to be responsible for the stability of the so-called "chiral" tubes, those which are neither armchair nor zigzag 11 .This thermodynamic model allowed to identify the most stable nanotube structure as a function of three parameters: the interface energies of the armchair or zigzag edge atoms and the temperature.Thus, distributions of the probability of nucleation of the tubes as a function of their chirality could be established, which compare favorably to the experimental results, but have a tendency to be too broad, probably because the growth kinetics were not taken into account.
Our goal here is thus to explore growth kinetics, using the simplest possible kinetic Monte Carlo (kMC) simulations on a fixed lattice unique to each helicity.A similar route has recently been proposed for studying the growth and etching of graphene flakes on catalytic substrates 12 .In our approach which is a development of the work by Dumlich and Reich 13 using more powerful computer simulation tools, we consider the growth from the point of view of the tube, with the underlying idea that the interplay between different edge geometries that are constrained by the helicity of the tube and the energies of the interface in contact with the catalyst, could lead to a chirality-selective growth.In this approach, the thermochemistry of the catalytic decomposition of the gaseous precursor, the diffusion of carbon atoms to the nanotube lip, and the amount of carbon dissolved in the catalytic nanoparticle (NP) that controls the growth mode 14 , are considered as elements of a global process, independent of the tube's helicity, that delivers carbon atoms at a given chemical potential to the tube edge and thus drives the growth or etching of the tube.
In the very simple thermodynamic model we proposed, the interface of a (n, m) tube was assumed to be perpendicular to the tube axis, with n + m carbon atoms in contact with the catalyst, 2m of them in armchair configuration and the remaining n − m in zigzag configuration.Such interfaces are indeed observed 11 , but recent High Resolution Transmission Electron Microscopy (HRTEM) observations of tube-catalyst interfaces after the synthesis show some variability of the interface structure.Nanotubes were synthesized using the floating catalyst CVD method 15 .
HRTEM samples were collected directly from the CVD reactor by placing the TEM grid (copper) with thin holey carbon supporting films (Agar Scientific Ltd.) on a filter for 20 seconds.A double aberration-corrected JEOL-2200FS microscope (JEOL Ltd.) was employed for high-resolution observation operated at an acceleration voltage of 200kV.In order to enhance the contrast of carbon nanotubes sitting on thin carbon films, HRTEM images were taken at slightly over-focussed conditions.A number of NPs were found isolated, either inactive or detached from the tube, but a few tens of tube-catalyst interfaces could be observed, showing that nanotubes may display more complex, oblique interfaces.Typical examples are displayed in Fig. 1.However, it should be noted that the structure of the interface might be modified during the cooling process, but, using a similar approach, Fiawoo et al. could identify different growth modes 16 that have been validated since then 14 .Such interfaces often correspond to the tube axis being off the center of the NP.This might help quantifying how oblique interfaces are.In the present example, a majority of the observed interfaces were oblique, but this probably depends on the CVD conditions, and a more detailed study should be performed.This points to an important observation, overlooked in Ref. 11: if the number of carbon atoms at the interface is kept equal to n + m, there is only one possible edge for zigzag tubes, while tubes with larger chiral angles θ can exhibit a large number of edges, including more oblique ones.
The most oblique interface of a (n, n) tube contains 2n − 2 zigzag and 2 armchair atoms.This asymmetry has consequences on both thermodynamic stability and growth mechanisms.This article is organised as follows.We start by defining our lattice model and the associated energies.We then present the kMC algorithm.Since its application in a grand canonical ensemble is not so common, we demonstrate the validity of our approach that makes it possible to simulate growth and etching of tubes.This enables us to calculate growth rates for a large set of 96 tubes with different chiralities, that depend on the five parameters of the energy model and the growth conditions (carbon chemical potential and temperature).We then try to disentangle the influence of each parameter and finally display three typical situations where a high selectivity based on large differences of growth rates can be expected.

II. THE MODEL
We model the CNT growth process using a lattice gas, where lattice sites are either empty or occupied by a carbon atom.Carbon atoms in the lattice have either three carbon neighbors when they are within the tube, forming "bulk" sites, denoted C 3 , or two carbon neighbors and one empty neighbor site, forming either zigzag C Z or armchair C A edge sites.The former have two neighboring C 3 sites, the latter have one C 3 and one two-fold coordinated sites.We also consider C atoms with only one C neighbor, attached to either C A or C Z edge sites, denoted C 1 A and C 1 Z , respectively.These high energy configurations approximate transition states for the addition and removal of atoms during the kMC simulations, as discussed below.Isolated C atoms with no C neighbor and other types of two-coordinated atoms, forming chains, are ruled out.Thus there are A , E 1 Z of these states are interface energies, as opposed to those of atoms at an open nanotube edge without contact to a catalyst.Such interface energies are very difficult to measure experimentally, but have been evaluated using DFT calculations for different catalysts.Most of these calculations were performed on small systems relaxed at zero temperature, and rely on assumptions regarding the structure of the catalyst NP and its surface.Catalysts probed in Ref. 17 (Au, Cu, Pd, Ni, Co, Fe) show a stronger interaction with zigzag tubes, but an investigation of a larger set of metals 18 indicates that this is not a general trend.Interestingly, using DFT-based molecular dynamics at high temperature for very small bimetallic clusters, Qiu and Ding 19 showed that the nanotube-metal interaction can be strong enough to drag the metal with the higher carbon affinity to the lip of tube, thus modifying the structure of the nanoparticle.Recent theoretical investigations 20,21 suggest that ordering or phase separation tendencies between armchair and zigzag species at the edge of the tube can be induced by the catalyst.We thus include an ordering energy term E A/Z = r ε A/Z , where r is the number of contacts between armchair and zigzag species and ε A/Z the energy per contact, to account for these tendencies.As illustrated in Fig. 2, the total energy of a configuration is the sum of two terms: a "bulk" term that depends on the energy and number of C 3 atoms and an interface energy term E Int that includes carbon atoms with less than three carbon neighbors.Setting the energy reference to zero for C 3 carbon atoms in the lattice, we end up with a model in which the total energy is reduced to the interface energy that writes: Note that all interface energy terms for 2-and 1-fold coordinated C atoms are positive, while a negative (positive) ε A/Z favors armchair-zigzag alternation (separation) at the edge.The typical range of E A and E Z is 0.0 − 0.5 eV/atom 11 , while calculations of the interface energies of various SWNTs in contact with the W-terminated (0 0 3) surface of a Co 7 W 6 alloy 20,21 led to positive values of ε A/Z , that would correspond to about +0.30 eV in our model.We note that none of these calculations dealt with the real conditions encountered during a CVD synthesis, in particular concerning the NP structural stability 22 and the role of carbon dissolution in the catalyst 23 .In the present approach, we consider these interface energies as parameters, allowing them a large flexibility to probe the limits of the model.
During the course of the kMC simulations, carbon atoms are incorporated in or removed from the lattice.The corresponding statistical mechanical framework is the so-called semi-grand canonical ensemble (sGC-kMC): the total number of occupied lattice sites fluctuates in a manner controlled by the temperature T and a chemical potential difference ∆µ C between empty and occupied sites in the lattice.Defining µ tube

C
as the chemical potential of a carbon atom incorporated in the tube and µ cat C as that of a C atom in or on the catalyst, defined by the thermochemistry of the decomposition of the precursor, we simulate the CVD process by identifying ∆µ C as the difference The growth and etching mechanisms are driven by ∆µ C .Equilibrium corresponds to ∆µ C = 0.
When ∆µ C = 0, the system is driven out of equilibrium with unbounded growth for ∆µ C > 0 and etching for ∆µ C < 0. We note that the free energy landscape of our simple model is closely related to the discrete state Brownian motor, a classic model of non-equilibrium statistical mechanics 24 .
empty sites interface with armchair and zigzag atoms grown tube Direction of growth

III. KINETIC MONTE CARLO ALGORITHM
We wish to understand how the growth rate for the chirality (n, m) is influenced by the energy parameters defined above and the CVD conditions (T , ∆µ C ) by sampling the structure of the tube/catalyst interface in Monte Carlo runs.Furthermore, our goal is to identify conditions of high selectivity for certain chiralities.We thus have to extend the classical kMC algorithm 25 to grand canonical conditions.This has already been done 26,27 in the context of fluid systems.Barriers correspond to rare events, that are sampled efficiently in kMC, giving thus access to the physically important phenomena.In the following, performing a set of kMC simulations under the same conditions, we compute growth rates of 96 SWNTs with chiralities ranging from (5, 4) to (18, 1) (diameter range: 6.0-14.8Å).KMC approaches have already been used in the context of SWNT growth, in particular in an early work by Maiti et al. 28 , and more recently by Li et al. 29 .The main differences in the present work are that we use very simple lattice systems to enable extremely fast simulations, which in addition allow the addition and removal of carbon atoms to exhibit both growth and etching.
A key challenge when building kMC models is the identification of all possible atomic level mechanisms, in general a highly non-trivial task 30 .In the present case the actual atomic scale growth process is extremely complex, dependent on the nature of the carbon precursor, the catalyst, and the specific CVD parameters.Using Density Functional Theory (DFT) calculations, Shibuta et al. 31 showed that the decomposition of a simple methane molecule on a Ni surface involves multiple steps and a subsurface burying of the dissociated carbon atoms.Decomposition of the different bonds of the more complex ethanol molecule has been shown to depend on the surface sites of either Fe, Co or FeCo catalysts 32 .Closer to the present context, Penev et al. 20 showed how carbon dimers can be attached at an armchair-zigzag contact, forming a pentagon/heptagon pair which can lead to a tube chirality change if it persists.
Such complex mechanisms cannot be taken into account in our fixed lattice simulations; we retain only a phenomenological free energy barrier for the addition or removal of two carbon atoms at selected edge sites.In the same spirit as the algorithm described by Voter 25 for the canonical ensemble, during our sGC-kMC simulations, moves are randomly selected from a set of insertion and removal events, transitioning from a basin state a to a basin state b.Carbon atoms can be deposed on any C A or C Z site that allows to complete a six-membered carbon ring upon addition of a second atom, as illustrated in Fig. 3.As the added carbon atom is a high energy, weakly metastable state, we assume that in a period of time much shorter than the typical residence time in any basin state, this atom either detaches or is stabilized through the addition of a second carbon atom.Whilst this is consistent with the degree of approximation inherent in our simplified model, we note that a refined or even completely distinct mechanism would only change the value of the free energy barrier for this process, not the catalogue of available moves in our model.As the transition state contains one additional atom, we assume that the free energy of this state changes by −∆µ C , which is expected to be valid provided that the model transition state has a free energy greater than the initial and final states.
At each kMC step an exhaustive list of such insertion and removal events is established and the energy of each of the corresponding barrier states is evaluated.These rules imply that for the basin states all n + m edge atoms are either of type C A or C Z .The edge configurations of the barrier states contain n + m − 1 atoms of these types and one mono-coordinated atom of type C 1 A or C 1 Z .Each full Monte Carlo step, i.e. from one basin state to the next, fills or empties two lattice sites, thereby creating or eliminating two C 3 atoms, respectively.We identify 14 possible barrier states, 5 for insertion and 9 for removal, as shown in table I.A graphical representation of these barrier states can be found in Fig. 10.Depending on nanotube chirality and model parameters not all of them necessarily occur in the simulations.In particular, in the case of m = 0, there are no barrier states, and therefore these tubes cannot grow within this model.In the case of m = 1, only a subset of 4 barriers are available, due to geometric constraints.The energy of the barrier states is calculated from table I and we notice that the energy barriers for carbon insertion and removal are in principle asymmetric.We neglect entropic contributions from thermal vibrations, meaning all free energies are determined solely by the energetic model and the chemical potential difference ∆µ C .Growth moves (C insertion) are associated with 5 barriers, with 3 independent parameters only, that are: ) and ε A/Z .Etching moves (C removal) correspond to 9 barriers, with 5 independent parameters only, that are the same as above plus We use the rejection-free algorithm as detailed by Voter 25 , deriving from the "n-fold way" algorithm by Bortz et al. 33 : From the full list of available barriers at each given MC step, a barrier i is randomly selected with a weight proportional to the corresponding rate constant, calculated as: where ∆E i is the energy barrier height, k B the Boltzmann constant and ∆N = +1 when adding or −1 when removing a C atom.The exponential term is simply proportional to the ratio of the thermodynamic probabilities of the barrier and the basin a states, in the semi grand canonical ensemble.The important factor that controls the process is In the following we refer to it as the "effective barrier".In principle, the attempt frequency ν 0 can be estimated from experimental data or molecular dynamics simulations.Since the purpose of this paper is to compare growth rates of CNTs with different chiralities, which will be unaffected by the absolute value of ν 0 , we take τ = 10 6 ν −1 0 as unit of time.At each step, the clock is advanced by drawing a random number from the distribution p (t) = k tot exp (−k tot t), with k tot = ∑ i k i .This means that, with this rejection-free algorithm, the selected move is systematically accepted, but the probability of selecting a particular barrier depends on its height.The corresponding time increment for a particular move, which defines the kinetics, depends on the height of all barriers.

A. Observables
Growth rates can be defined either as the mass uptake, or as the tube length added, both scaled per unit of time τ.When dividing the latter by the tube diameter, the two are proportional.The former therefore favors smaller diameter tubes, and the latter larger diameter ones.Except for tubes with m = 0, which do not grow, and m = 1 with four available barriers, independent of diameter, the number of available barriers is roughly proportional to the circumference of the tubes.The clock advances proportional to k −1 tot , and therefore roughly inversely proportional to the tubes circumference for m > 1.Therefore, for these tubes, the length of tube added per unit of time is approximately independent of the diameter.As a consequence, significant differences in growth rate with respect to diameter cannot be expected and our calculated growth rates mostly depend on chiral angles.In this model, all tubes with m ≥ 1 grow to some extent.The chiral angle of a (n, m) tube is defined by θ = arctan √ 3m m+2n .Values range between 0 • for zigzag (n, 0) and 30 • for armchair (n, n) tubes.We group the chiral angles of the tubes in 2 • bins, and average the corresponding growth rates to analyze the selectivity.We define the selectivity S(θ ) as the growth rate in terms of tube length G(θ ) divided by the sum of all growth rates of the set of all 15 bins.
IV. RESULTS

A. Validation of the modeling
Before analyzing the effect of the model parameters on growth rates, we establish a baseline result which can be compared to an analytical expression for the growth rate.For this, two-fold coordinated atoms at the interface are assigned E A = E Z = 0.15 eV, one-fold coordinated atoms 8 eV, and the ordering energy ε A/Z = 0.0 eV.In this case, all energy barriers are equal, ∆E i = 0.65 eV.The temperature is set to T = 1000 K and the chemical potential difference to ∆µ C = 0.1 eV.In the following we use these values in simulations of 10 5 steps, where not indicated otherwise.
With these baseline parameters, tubes with large chiral angles are favored, with a relatively low selectivity (see Fig. 4).As expected, tube diameter does not affect growth rates very much.
More interestingly, because all barriers and hence all time increments are drawn from the same probability distribution (assuming equal n + m), the faster growth on the armchair side results from the larger number of available of edge sites allowing the formation of a new hexagon by adding two carbon atoms.Ultimately, because each dimer addition at the edge leads to exchanges of armchair and zigzag edge atoms, this is related to the larger configurational entropy of the edges on the armchair side.We then analyzed the effect of the chemical potential on growth rate.Except for zigzag tubes, which are frozen at any value of ∆µ C , all tubes are decaying for ∆µ C < 0, and growing for ∆µ C > 0.
At ∆µ C = 0, the tube height fluctuates around its initial position, with a constant average value.
In the baseline case, the fastest growing tubes are also those that are etched away the fastest, and the ratios of the growth rates for the different chiralities remain constant as the chemical potential is changed.Figure 5(a) shows the average growth rate of the set of 96 CNTs as a function of the chemical potential.For this baseline parameter set it can be shown that (see appendix B), where we find N = 5.207±0.004, the average number of insertion and removal sites.With varying ∆µ C , we confirm that the SWNT height obeys the fluctuation-dissipation theorem as ∆µ C → 0, see appendix C, demonstrating that our simulations are sampling the correct thermodynamic ensemble.
In order to evaluate the effect of temperature on the growth rate, simulations in the temperature range of 500-1500 K were carried out in the case of the baseline parameters with ∆µ C = 0.1 eV.
As with the dependence on ∆µ C , the temperature did not affect the ratios of the growth rates of the different chiralities.In this series of simulations, where baseline parameters have been used, neither temperature nor chemical potential significantly affects the proportionality constants between the growth rates of the different chiralities.However, looking at the components of the energy barriers in table I and the expression of the rate constant [eq.( 1)], we can see that other sets of parameters may exist where the growth rates are not always proportional.In particular, when energy barriers are not equal, growth and etching rates can be asymmetric with respect to chemical potential.

B. Sensitivity of the model parameters
In the following, we try to understand the influence of the interface energy parameters on the (n, m) dependent growth rates, and to identify conditions that could lead to a chirality-selective growth of SWNTs.As discussed in section III, these energies set the heights of the barriers that can be encountered during the simulations.
We begin the discussion with the influence of the interface energy E Z of zigzag atoms, keeping the other energy parameters at their baseline values.The evolution of the growth rates as a function of chirality and E Z can be appreciated from the animated Fig. 4. Results are quantified in Fig. 6 that shows the chiral angle of the fastest growing tubes as a function of E Z and the corresponding selectivity.At vanishing E Z , the growth rate of (n, 1) tubes is the largest and the growth is highly selective, because (n, 1) tubes can grow with N Z = n + m − 2 = n − 1 zigzag edge atoms that are energetically favored, and only N A = 2 armchair edge atoms.Moreover, growth is mostly limited to these (n, 1) tubes.As can be seen from Fig. 7(a), the growth proceeds by spiraling around the tube axis which is the only growth process available for tubes with m = 1.Only insertion barriers 1 and 4, and removal barriers 13 and 14 are utilized.All four barriers are available at every kMC step.However, with the present choice of parameters, barrier 4 has an effective height of only 0.4 eV, compared to the second lowest barrier 13 with 0.6 eV.Therefore barrier 4 is chosen in most of the MC steps.The time increment is small as it is dominated by this low energy barrier, which in turn leads to high growth rates.For other chiralities, the growth is also dominated by barrier 4, which is however, not always available.In these cases, more complex growth patterns emerge that are efficiently sampled in kMC.However, the available barriers are higher and time advances by a large increment, which significantly decreases overall growth rates for any tube with More available barriers at higher chiral angles directly translates into a smaller time increment which leads in turn to higher growth rates for these chiralities.Increasing E Z further, gradually decreases the growth rate for tubes with large chiral angles, so that the highest growth rates occur for tubes with chirality close to (2n, n), however at a rather low selectivity.Having discussed the influence of the energy of zigzag atoms, we now turn to the energy of other types of interface atoms.As barrier states always have one one-fold coordinated atom and n + m − 1 two-fold coordinated atoms, all barrier energies depend on differences involving 1-and 2-coordinated atoms.Therefore, the effect of the energy of 1-coordinated atoms E 1 X is opposite to that of 2-coordinated atoms E X , where X stands for A or Z.In particular, a small E 1 Z leads to a similar growth rate profile as obtained with a large E Z and vice versa.Similarly, the effect of armchair and zigzag energies on growth rates are opposite.Again, this can easily be understood from the barriers.What affects the growth rates is the relative energy of the different barriers.In particular, barriers 1, 3, 7, 9, 12, and 14 decrease with E Z , while barriers 2, 4, 5, 6, 8, 10, 11, and 13 decrease with E A .This means both E A and E Z change the relative height of the two groups of barriers, but in opposite directions.The types of growth rate maps obtained while varying these other interface energies (E 1 Z , E A , and E 1 A ) are not qualitatively different from those obtained while varying E Z .

At intermediate values of zigzag energy (E
Finally, we consider the effect of the ordering energy parameter ε A/Z on the growth rate profiles.Figure 8 shows the chiral angle of the fastest growing tube (left y-axis) and the selectivity associated with this growth (right y-axis) as a function of ε A/Z for two different values of E A .
For negative values of ε A/Z , growth profiles are qualitatively very similar to the baseline case, i.e. armchair tubes have the fastest growth, but with relatively low selectivity.For positive ε A/Z , and E A = 0.15 eV, (m, 1) tubes have again the fastest growth, but for E A = 0.05 eV near armchair (n, n − 1) tubes grow the fastest with a higher selectivity.We are again dealing with a type of spiral growth, but dissimilar from the one encountered for the (n, 1) tubes, see Fig. 7(b).This growth process also occurs in a somewhat similar fashion for armchair tubes, however at each added new layer, a new kink needs to which slows down the growth process.The barrier associated with this growth mode is barrier 3, with an effective height of 0.25 eV, which is much lower than the next barrier 4 at 0.45 eV.Other tubes grow much slower, because barrier 3 is not always available, which makes the time increments much larger.

C. Examples of high selectivity
Having understood how the parameters of the model affect growth rates, we explore the possibility to achieve a high selectivity on the chiral angle distributions.Figure 9 shows growth rate profiles for three cases of chiral angle selective growth.As before, all the following kMC runs are performed at 1000 K, however, the energy parameters are given a larger flexibility.
Panel (a) shows that a faster growth of (n, 1) tubes is obtained for E A = 0.172 eV, E Z = 0.05 eV, E 1 A = 0.8 eV, E 1 Z = 1.5 eV, ε A/Z = 0.1 eV, and ∆µ C = 0.1 eV.There is only one possible edge structure for (n, 1) tubes, with (n + m − 2) zigzag and two armchair edge atoms, implying that spiral growth is their only possible mechanism.Since (n, 0) tubes do not grow, (n, 1) tubes are competing with tubes with larger chiral angles, which offer an increasingly large number of possible edges, by transforming armchair into zigzag edges and tilting the interface with respect to the tube axis.Because zigzag edges are more stable with these parameters, all (n, m) tubes with m > 1 tend to maximize their number, while mostly (60% of the configurations) keeping at least two armchair pairs instead of only one.The reason is probably the entropy gain that stabilizes the interface, while the two partial spiral developments hinder each other.If an anti-armchair motif (two neighboring zigzag edges with different inclinations) is present, which is often the case for m > 1, an armchair pair can be added with barrier 5 and subsequently be removed using barrier 12.This leads to oscillations that pause growth over extended periods of time.The mechanism occurs only if barrier 4 (effective height 0.21 eV) is not available, because these barriers have an effective height of 0.41 eV and 0.85 eV which is much higher.These high barriers also make the time increment much larger, which slows growth even more.These are the reasons why all (n, m) tubes with m > 1 have similar growth rates about one thousand times smaller than (m, 1) in the present conditions.We note however that E Z is generally larger than E A on metal catalysts 34,35 , explaining why this calculated preference (n, 1) tubes is rarely obtained.

V. DISCUSSION
A central assumption in this study is to consider ideally crystalline tube structures, without defects.Once created during the growth process 38 , defects can be healed in the vicinity of the catalyst 39 or remain present.Whether defect healing is correlated with the tube structure 40,41 is still an open question, and we assume that the presence of defects has no influence, other than detrimental, on the mechanisms fixing the chirality during growth.Since our goal is to answer the question if growth kinetics can lead to a significant chiral selectivity, we have to consider ideal crystalline tube structures and compare their growth rates: ideal tube structures are formed with carbon hexagons only, and we consequently consider only Monte Carlo moves that add or remove carbon dimers at the edge of the tube.A consequence is the somewhat oversimplified definition of the barrier states.On the upside, this avoids discussing detailed atomistic mechanisms that remain inaccessible experimentally, and the energies of one-fold coordinated species at the edge E 1 Z and E 1 A can be seen as variables controlling barrier heights that are necessarily present in the process, though essentially unknown, but of prime importance on the kinetics.
Our simulations show that SWNT growth is determined by complex processes involving the relative stability of armchair and zigzag edges, the degrees of freedom of the edge structure, and the associated configurational entropy, as well as the energy barriers controlling the incorporation of incoming carbon atoms.The entanglement of these contributions makes it difficult to develop a rational approach to experimentally achieve a high selectivity, but we expect this contribution to promote it.The spiral growth process envisioned by Ding et al. 34 is indeed sometimes observed.However, the conclusions were drawn too hastily.In fact, not only (n, n − 1) but also (n, 1) series may grow fast through this process, but under some conditions more disordered growth, involving a larger set of carbon incorporation barriers proves more efficient than pure spiral growth.
In our previous modeling 11 , we assumed that the interface was perpendicular to the tube axis, with fixed numbers of armchair and zigzag edge atoms.Since oblique interfaces have been observed experimentally, the natural step forward is to assume that edges simply have n + m edge atoms.Allowing for even longer edges would lead to unphysically large interface meandering.
In the present model, zigzag tubes cannot nucleate and grow.DFT calculations by Ding et al. 34 show that the probability of nucleating zigzag tubes is very low.Another consequence is that the configurational entropy of the edge is larger on the armchair side.A full development of the thermodynamic aspects of these assumptions should be performed to obtain new (n, m)-dependent probability distributions, that, combined with the growth rates calculated here would produce narrower diameter dependent distributions.Maps showing the abundance of each tube chirality could then be more significantly compared with experimental data than the purely kinetic data presented here.Such calculations are beyond the scope of this introductory paper.

VI. CONCLUSIONS
We have developed a simple yet powerful kinetic Monte Carlo modeling of SWNT growth that allows comparison of growth rates of large sets of tubes with various chiralities.Significant differences in growth rates can be expected by tuning the interfacial energies.This growth selectivity of kinetic origin complements the thermodynamic selectivity 11 operating at the nucleation stage.In principle, interfacial energies depend on the growth mode 14 and the adhesion energy of the tube to the catalyst, which can be evaluated by DFT calculations [17][18][19]42 . In ractice, the interfacial properties should be adjusted by an adequate choice of catalyst and experimental conditions.Direct validation of our results is difficult due to the scarcity of measurements of growth rates of individual tubes 43,44 .However, the progress currently made in this direction 45 is a strong motivation for this work.
Beyond the general analysis of the growth kinetics which critically depends on the identification of all possible barriers, we have identified growth parameters (interfacial energies, carbon chemical potential, temperature) leading to experimentally observed selectivity patterns.Baseline parameters for which all insertion and all etching barriers are equal lead to a slight selectivity towards large chiral angles, commonly observed for tubes grown with Fe in perpendicular mode 14 .Sharper near armchair selectivity is obtained for interfacial energy parameters E A < E Z and when energy barriers 3 and 7 are the weakest.The first condition corresponds to (E adh A + E db A ) < (E adh Z + E db Z ), which is observed for many catalysts 17 .In such a case, the difference E adh A − E adh Z between the adhesion energies of armchair and zigzag edges is less than the difference in dangling bond energies E db Z − E db A ≈ 1.1eV .Finally, we have identified conditions for growing tubes with a chiral angle around 19 • , which come at the price of a low yield, as observed experimentally, because of a rather inefficient growth mechanism.
Our simulations also highlight the importance of including both insertion and removal mechanisms in the kMC algorithm.With the non-selective baseline parameters, growth and etching are symmetrical for all tube chiralities, as observed experimentally 44 .However, our calculations  green symbols) and the mobility v(∆µ C )/(∆µ C /l) in the limit ∆µ C → 0 (small, purple symbols), which are predicted to be equal by the fluctuation-dissipation theorem, eq.(C1).Circles and squares correspond to kMC simulations using the "baseline" parameter set and E Z = 0.5 eV, respectively.

FIG. 1 .
FIG. 1. (Color online) Different carbon nanotube-catalyst interfaces: (a-d) HRTEM images of carbon nanotubes (yellow) attached to catalyst nanoparticles.The interfaces are highlighted in red.In some cases the center of the SWNT is not aligned with the catalyst which leads to oblique interfaces.(e) Armchair tubes, here (6, 6), can have an oblique rim of n + m = 12 atoms with up to 2n − 2 = 10 zigzag (green) and 2 armchair (pink) atoms.(f) An oblique interface of the zigzag (12, 0) tube has more than n + m = 12 atoms, here 16.

FIG. 2 .
FIG. 2. (Color online) Schema for the lattice kMC simulations of the growth of carbon nanotubes.Carbon atoms in zigzag and armchair configuration in contact with the catalyst (not represented) at the edge of the tube are shown as green and pink balls, respectively.Atoms in the bulk part of the tube are represented as gray balls and the empty part of the lattice is drawn as white sticks.
FIG. 3. (Color online) Basins (initial and final state) and barrier states encountered in the kMC algorithm.The barrier states always contain a mono-coordinated atom (represented as a blue sphere) at the interface with the catalyst, while the basin states contain only zigzag and armchair atoms (green and violet spheres, respectively).Energy parameters are chosen in such a way that only positive effective barrier heights are encountered in the simulations.

FIG. 4 . 1 A = E 1 Z
FIG. 4. (Color online) Growth rates for individual chiralities in the "baseline" case : E A = E Z = 0.15 eV, E 1 A = E 1 Z = 0.8 eV, T = 1000 K and ∆µ C = 0.1 eV.Zigzag (n, 0) tubes do not grow, while growth rates increase with chiral angle and only very weakly depend on the diameter.The animation shows maps for 0.0 eV ≤ E Z ≤0.5 eV including the corresponding barrier heights (Multimedia view).
Z ≈ 0.15 eV), we recover the baseline condition, where armchair tubes grow the fastest, but with a relatively low selectivity (about a third of what was achieved with the tubes with m = 1).The transition between favoring tubes with m = 1 to n = m is abrupt at E Z ≈ 0.05 eV.For E Z ≈ 0.15 eV, only two different values of effective barrier energies exist: all insertion barriers are at 0.55 eV and all removal barriers at 0.75 eV.The insertion barriers are therefore selected completely randomly.The number of available sites grows with chiral angle, as insertion cannot occur at two neighboring zigzag edges with the same inclination.

FIG. 6 .
FIG. 6. (Color online) Favored chiral angle (purple, left y-axis) and growth selectivity (green, right y-axis) as a function of the interface energy of zigzag atoms.

FIG. 7 .
FIG. 7. (Color online) (a) Spiral growth for a (11, 1) tube at vanishing E Z , (b) Spiral growth for a (7, 6) tube with E A = 0.05 eV and ε A/Z = 0.15 eV, and (c) Non-spiral growth of a (8, 4) tube with parameters corresponding to Fig. 9(c).Armchair and zigzag atoms are represented as violet and green spheres respectively.The animations show the growth mechanisms in the three situations, but it is important to note that the sequence of configurations does not include the time delays introduced by the kMC algorithm and hence do not give any hint on the kinetics (Multimedia view).

FIG. 8 .
FIG. 8. (Color online) Effect of ordering energy on favored chiral angle (left y-axis, purple symbols) and growth selectivity (right y-axis, green symbols).Full circles and open squares and correspond to E A = 0.05 eV and E A = 0.15 eV, respectively.
FIG. 10. (Color online) (a-n) Barrier states (1-14) encountered in the kMC simulations.The barrier states always contain a mono-coordinated atom (either attached to an armchair atom [blue sphere], or attached to a zigzag atom [orange sphere]) at the interface with the catalyst.Zigzag and armchair atoms are shown as green and violet spheres, respectively.

FIG. 11
FIG. 11. (Color online) The diffusivity divided by the thermal energy, D/k B T , obtained at ∆µ C = 0 (large,

TABLE I .
List of energy barriers occurring in the kMC simulations given as a function of the difference in the number of geometric features of the interface.