## Abstract

Experimental results (Huang *et al*.) indicated that nanotwinned diamond (nt-diamond) has unprecedented hardness, whose physical mechanism has remained elusive. In this report, we categorize interaction modes between dislocations and twin planes in nt-diamond and calculate the associated reaction heat, activation energies, and barrier strength using molecular dynamics. On the basis of the Sachs model, twin thickness dependence of nt-diamond hardness is evaluated, which is in good agreement with the experimental data. We show that two factors contribute to the unusually high hardness of nt-diamond: high lattice frictional stress by the nature of carbon bonding in diamond and high athermal stress due to the Hall-Petch effect. Both factors stem from the low activation volumes and high activation energy for dislocation nucleation and propagation in diamond twin planes. This work provides new insights into hardening mechanisms in nt-diamond and will be helpful for developing new superhard materials in the future.

## INTRODUCTION

Diamond, the hardest material available (*1*), is of great interest both scientifically, for understanding the origin of its outstanding mechanical properties (*2*), and technologically, for its direct industrial and research applications as cutting and polishing tools, diamond knives, diamond anvil cells, etc. (*3*). Therefore, numerous experimental and theoretical studies have been performed to understand its origin of superhardness (*4*, *5*) and improve its hardness further (*6*–*9*). Similar to metals and alloys, diamond can be hardened by grain refinement. The hardness of synthetic nanograined diamond (110 to 140 GPa) is higher than the average hardness of single crystal diamond (96 GPa) (*7*). Recently, nanotwinned diamond (nt-diamond) synthesized by compressing onion carbon has achieved hardness values more than twice that of natural diamond (180 to 200 GPa) (*8*, *9*), setting a new world record. Can the hardness be further increased? Understanding the origin of hardness in nt-diamond is an important scientific problem, as it may provide a new strategy for designing new materials with even higher hardness.

Although it has been speculated that the unprecedented hardness in nt-diamond may be due to combined Hall-Petch effect (*10*, *11*) and quantum confinement effect (*12*–*14*), alternative mechanisms have also been proposed (*15*, *16*). Diamond is a covalent material, and its elastic and plastic properties are largely different from those of metals. As carriers of plastic deformation (*17*), dislocations in diamond must also be different from those in metals. The Hall-Petch effect is an experimental rule deduced from metals (*18*). Is it valid to apply this rule to covalent materials? To answer these questions, a detailed study on interactions between dislocations and twin planes in nt-diamond is necessary. Here, we investigate these interactions using molecular dynamics (MD) simulations. Our calculated results indicated that two factors contribute to the unprecedented hardness in nt-diamond: high lattice frictional stress by nature of diamond and high athermal stress due to the Hall-Petch effect. This work provides a new insight into the nt-diamond hardening mechanism, which will be helpful for studying superhard materials in experiments.

## RESULTS

### Dislocation types and slip modes in nt-diamond

The structure of diamond can be envisioned as two interpenetrating face-centered cubic (fcc) lattices. Dislocations in diamond lie primarily along the < 110 > directions (*19*), which are associated with deep energy troughs, requiring large Peierls stresses to be mobile. There are two inequivalent {111} slip planes: the shuffle-set and the glide-set (*20*, *21*). Accordingly, there are four groups of perfect dislocations, namely, glide-set 0° (screw) dislocations, glide-set 60° dislocations, shuffle-set 0° dislocations, and shuffle-set 60° dislocations. Similar to fcc metals, a glide-set perfect dislocation can disassociate into two glide-set partial dislocations. Usually, a glide-set 0° perfect dislocation can disassociate into two glide-set 30° partial dislocations, and a glide-set 60° perfect dislocation can disassociate into a glide-set 30° partial dislocation plus a glide-set 90° partial dislocation. In contrast, for shuffle-set dislocations slipping along the <112> direction, there is no stable stacking fault (SF) and the corresponding unstable SF energy is very high. Therefore, shuffle-set perfect dislocations cannot disassociate into partial dislocations (*22*). This results in six types of dislocations in diamond: glide-set 0° perfect dislocations, glide-set 30° partial dislocations, glide-set 60° perfect dislocations, glide-set 90° partial dislocations, shuffle-set 0° perfect dislocations, and shuffle-set 60° perfect dislocations.

As illustrated in Fig. 1A [simplified on the basis of the micrographs reported in (*8*)], nt-diamond is composed of near-parallel nanoscale twin lamella in submicrometer-sized grains. Across a given twin boundary (TB), the slip systems can be expressed as a combination of two opposite Thompson tetrahedrons (Fig. 1B). Tetrahedrons ABCD and A^{T}B^{T}C^{T}D^{T} represent the Burgers vectors of dislocations in the parent crystal and its twined counterpart, respectively. In either case, Burgers vectors of perfect dislocation are parallel to the edges of the tetrahedron (such as AB, AC, BC, etc.). The partial dislocations are denoted by Roman-Greek pairs (such as Aδ, Cδ, Bδ, etc.). According to dislocation lines and their kink orientations, we divide dislocation slip modes into three types in nt-diamond, all of which are plotted in Fig. 1B (*23*). In mode I slip (slip transfer mode), dislocation lines are parallel to the twin planes, and the kink direction is nonparallel to the twin planes. In this mode, the dislocations can penetrate TBs into the adjoining twin lamellae or be absorbed into the twin plane. In mode II slip (confined layer slip mode), both dislocation lines and the kink directions are nonparallel to the twin planes. Dislocation motion is confined inside the independent twin lamellae. In mode III slip (paralleled to twin plane slip mode), both dislocation lines and kink directions are parallel to the twin planes. Both slip planes and slip directions are parallel to twin planes.

Critical resolved shear stresses (CRSSs) for the three slip modes are investigated. For slip mode I, when a dislocation reaches a twin plane, it may react with the twin plane and produce a new dislocation. The energy required for the production of the new dislocation results in an additional dislocation motion resistance, which, in essence, is the origin of the Hall-Petch effect. The CRSS for mode I slip therefore can be deduced by calculating the threshold stress of the dislocation reaction. For slip modes II and III, dislocation motion processes are confined within either the twin planes or individual grains, respectively, and their CRSSs can be deduced from the well-known dislocation motion theory.

### Reaction heat, activation energy, and barrier strength for the slip transfer mode

Glide-set 0° perfect dislocations. All four possible interactions between a glide-set 0° perfect dislocation and twin planes by slip mode I have been investigated. The interaction schematics and the corresponding reaction heat and activation energy are plotted in Fig. 2 (A and B). The interactions involve both perfect and dissociated dislocations (Fig. 2A). For the perfect dislocation interacting with a twin plane, there are two reaction pathways: BC → Bδ + δC and BC → B^{T}α^{T} + α^{T}C^{T}. When the perfect dislocation BC in the BCD plane reaches the twin plane ABC, BC dissociates into two partial dislocations (Bδ and δC) in ABC. At the same time, a twin plane migrates in the opposite direction of dislocation motion. By comparing the system energies before and after the dissociation, the corresponding reaction heat can be calculated. In addition, activation energy of the reaction has been calculated by building a dislocation kink motion path as shown in fig. S2. On the basis of the calculated kink energy and activation energy of kink migration (table S1), the obtained reaction heat and activation energy are −0.9 eV/Å and 6.2 eV, respectively. The activation energies for the dislocation reactions vary with applied shear stress (fig. S3). We define the barrier strength as the threshold stress of dislocation reaction with TB when the activation energy is equal to zero (*24*). From the shear stress–dependent activation energies thus obtained, the barrier strength is calculated to be 43.7 GPa. All the values are listed in Table 1 and plotted in Fig. 2B. Similarly, for the dislocation reaction BC → B^{T}α^{T} + α^{T}C^{T}, the glide-set 0° perfect dislocation dissociates into two glide-set 30° partial dislocations (B^{T}α^{T} and α^{T}C^{T}) separated by an SF. The calculated reaction heat, activation energy, and barrier strength are −0.9 eV/Å, 10.2 eV, and 31.7 GPa, respectively.

For a dissociated perfect dislocation interacting with a twin plane, there are two reaction pathways: Bα + αC → Bδ + δC and Bα + αC → B^{T}α^{T} + α^{T}C^{T}. For the former reaction, the partial dislocation Bα slips along the twin plane and becomes Bδ, leaving with a stair-rod dislocation (δα) in the twin plane. The trailing 30° partial dislocation αC then slips after the leading dislocation, capturing the stair-rod dislocation δα to form a glide-set 30° partial dislocation δC which slips in the twin plane. The calculated reaction heat, activation energy, and barrier strength are 0 eV/Å, 6.9 eV, and 53.7 GPa, respectively. For the latter reaction, the leading dislocation Bα penetrates the twin plane and becomes B^{T}α^{T}, leaving a stair-rod dislocation (α^{T}α) in the twin plane. The trailing 30° partial dislocation αC then slips after the leading dislocation, capturing the stair-rod dislocation α^{T}α and forming a glide-set 30° partial dislocation α^{T}C^{T} which slips into the adjacent twin domain. The calculated reaction heat, activation energy, and barrier strength are 0 eV/Å, 6.7 eV, and 52.1 GPa, respectively. We use the barrier strength as the criterion to identify the most possible dislocation interaction. The lower the barrier strength, the more possible for the dislocation interaction to occur. On the basis of this criterion, perfect glide-set screw dislocation BC interacting with twin planes is the easiest mode, as its barrier strength is the lowest among the four dislocation interaction modes.

Glide-set 30° partial dislocations. For the glide-set 30° partial dislocation Bα interacting with a twin plane (Fig. 2A), two reaction pathways exist: Bα → Bδ + δα and Bα → B^{T}α^{T} + α^{T}α. For Bα → Bδ + δα, the 30° partial dislocation Bα slips in the twin plane and becomes Bδ, causing the TB to migrate in the opposite direction of dislocation motion and leaving a stair-rod dislocation δα (*b* = 1/6<110>) in the twin plane. The calculated reaction heat, activation energy, and barrier strength (Table 1 and Fig. 2B) are 1.5 eV/Å, 6.9 eV, and 53.7 GPa, respectively. For the reaction Bα → B^{T}α^{T} + α^{T}α, the 30° partial dislocation Bα penetrates the twin plane and slips into the adjacent twin domain to become B^{T}α^{T}, leaving a stair-rod dislocation α^{T}α (*b* = 2/9<111>) in the twin plane (*25*). The reaction heat, activation energy, and barrier strength for this reaction are 2.4 eV/Å, 6.7 eV, and 52.1 GPa, respectively. Because the magnitude of the Burgers vector for α^{T}α is greater than that of δα, the reaction heat is higher than that of Bα → Bδ + δα. However, the energy barriers and barrier strength are almost the same for these two kinds of reactions.

Glide-set 60° perfect dislocations. The reaction of glide-set 60° perfect dislocation with a twin plane can occur in two ways: BD → Bδ + δD and αD + Bα → 2α^{T}δ + D^{T}α^{T} + B^{T}α^{T} (Fig. 2C). For BD → Bδ + δD, when a glide-set 60° perfect dislocation BD reaches the twin plane, the interaction releases a partial dislocation Bδ in the twin plane, leaving an additional Frank partial δD, and at the same time the TB migrates in the opposite direction of Bδ motion. The reaction heat, activation energy, and barrier strength (Table 1 and Fig. 2D) are 2.7 eV/Å, 9.9 eV, and 55.3 GPa, respectively. These values are higher than those of the reaction of a glide-set 0° perfect dislocation with the twin plane. For αD + Bα → 2α^{T}δ + D^{T}α^{T} + B^{T}α^{T}, a glide-set 60° perfect dislocation BD first dissociates into a glide-set 30° partial dislocation Bα and a glide-set 90° partial dislocation αD, separated by an SF. The leading 30° partial dislocation Bα penetrates the twin plane by the reaction Bα → B^{T}α^{T} + α^{T}α; the trailing 90° partial dislocation αD reacts with the stair-rod dislocation α^{T}α by the dislocation reaction αD + α^{T}α → 2α^{T}δ + D^{T}α^{T}. Therefore, the total reaction can be expressed as αD + Bα → 2α^{T}δ + D^{T}α^{T} + B^{T}α^{T}. The corresponding reaction heat, activation energy, and barrier strength (Table 1 and Fig. 2D) are 3.0 eV/Å, 6.7 eV, and 52.1 GPa, respectively. Because two stair-rod dislocations are produced in this reaction, the heat of the reaction is higher than that of the glide-set 0° perfect dislocation, as described in the “Glide-set 0° perfect dislocations” section. The energy barrier in this case is identical to that described in the “Glide-set 0° perfect dislocations” section.

Glide-set 90° partial dislocations. Two possible interactions for a glide-set 90° perfect dislocation with a twin plane by slip mode I are αD → Aδ and αD → αδ + α^{T}δ + D^{T}α^{T} (Fig. 2C). For reaction αD → Aδ, when a glide-set 90° partial dislocation αD reaches the twin plane, it becomes equivalent to a glide-set 90° partial dislocation Aδ in ABC, causing twin plane migration. Since both αD and Aδ are the same type of partial dislocations, the net total reaction heat is 0. The activation energy and barrier strength are calculated to be 3.9 eV and 37.2 GPa, respectively. For the reaction αD → αδ + α^{T}δ + D^{T}α^{T}, first, a 90° partial dislocation αD dissociates into a stair-rod dislocation αδ and a Frank partial δD. The Burgers vector of δD is identical to that of D^{T}δ, which can further dissociate into a 90° partial dislocation in the adjacent twin, along with a stair-rod dislocation in the twin plane: D^{T}δ → D^{T}α^{T} + α^{T}δ. The corresponding reaction heat, activation energy, and barrier strength (Table 1 and Fig. 2D) are 2.9 eV/Å, 4.6 eV, and 49.3 GPa, respectively. Because of two stair-rod dislocations produced in this reaction, its reaction heat is almost twice that of the 30° partial dislocation. The activation energy is lower than that of the 30° partial dislocation interaction described in the “Glide-set 30° partial dislocations” section.

Shuffle-set 0° perfect dislocation. The shuffle-set 0° perfect dislocation BC is always maintained as a full dislocation with pure screw properties. Therefore, the reaction of the shuffle-set 0° perfect dislocation with the TB (in Fig. 2E) has only two pathways: cross-slip reaction BC → BC and penetration through the twin plane into the adjacent twin domain BC → B^{T}C^{T}. The reaction heat for both reactions is 0 eV/Å. As shown in Fig. 2F, the activation energy and barrier strength (Table 1 and Fig. 2F) for BC → BC are 7.8 eV and 24.1 GPa, respectively, and those for the BC → B^{T}C^{T} reaction are 8.2 eV and 19.2 GPa, respectively. Thus, for both reaction pathways, the energy and shear stress requirement is the same.

Shuffle-set 60° perfect dislocation. The reactions of a shuffle-set 60° perfect dislocation with a twin plane (Fig. 2E) are BD → BC + CD and BD → B^{T}C^{T} + CD. For both reactions, when a shuffle-set 60° perfect dislocation BD reaches the twin plane, a new shuffle-set 60° perfect dislocation CD appears. Since Burges vectors produced for these two kinds of reactions are the same, their reaction heat and activation energy (Table 1 and Fig. 2F) are identical, that is, 7.1 eV/Å and 16.3 eV, respectively. The barrier strength of these two dislocation reactions is also the same (47.7 GPa).

The net Burgers vector must remain unchanged after a dislocation reacting with a twin plane. For the aforementioned six types of dislocations, their reactions with a TB also follow this rule. For example, for both shuffle-set and glide-set screw perfect dislocations, their reactions with a twin plane are carried out by cross slip, and there is no new dislocation formed. However, for glide-set 60° perfect, glide-set 30° partial, glide-set 90° partial, and shuffle-set 60° perfect dislocations, additional dislocations must be produced when they react with TBs so that the net Burgers vector remains unchanged after the reaction (*25*).

## DISCUSSION

### CRSS for slip mode I

To obtain the CRSS (τ_{css}) for slip mode I, a dislocation pileup model is used, whose schematics are plotted in the inset of Fig. 3A. According to this model, the CRSS can be expressed as the following (*26*)(1)where τ_{0} is the lattice frictional stress, G is the shear modulus, *b* is the magnitude of the Burgers vector, λ is the twin thickness, and τ_{TB} is the barrier strength of the dislocations when reacting the TBs. Both modulus and stress are in gigapascal, and all length parameters are in nanometers. The shear modulus of diamond is 540 GPa. The barrier strength τ_{TB} of shuffle-set 0° perfect dislocation penetrating into the adjacent twin is lowest in all dislocation reactions (19.2 GPa as shown in Table 1 and table S2). On the other hand, by calculating the stress-dependent activation energy for the shuffle-set 0° perfect dislocation BC slipping in a perfect diamond crystal, lattice frictional stress τ_{0} is obtained to be 10.3 GPa. Our result is consistent with the CRSS (10.7 GPa) of perfect diamond calculated from the diamond hardness (96 GPa), according to the relation τ_{0} = H/9. All these parameters are summarized in table S3.

By using Eq. 1 and the parameters summarized in table S3, a twin thickness–dependent CRSS for slip mode I can be expressed as(2)As shown in Fig. 3A, the CRSS for slip mode I increases with decreasing twin thickness.

### CRSS of slip mode II

Slip mode II is a dislocation motion confined by two twin planes with a defined twin thickness, so its CRSS can be described by the confined layer slip mode (*27*), whose schematics are also plotted in the inset of Fig. 3A. Here, the CRSS is equal to the extra stress needed to overcome the increased dislocation energy (Δ*W*_{Dis}). Therefore, τ_{css} of slip mode II can be expressed as (*27*)(3)Here, θ is the angle between the slip plane and the twin plane, and λ is the twin thickness (in nanometers). The dislocation energy can be expressed as follows(4)where φ is the angle between the dislocation line and the Burgers vector, ν is the Poisson ratio, and α is the dislocation core parameter (*25*). Again, all length-related parameters in Eqs. 3 and 4 are in nanometers. On the basis of Eqs. 3 and 4, the CRSS of slip mode II for shuffle-set 0° perfect dislocation AD motion is(5)In addition, on the basis of the material parameters of diamond (table S3), the resulting CRSS is(6)The calculated CRSS as a function of twin thickness λ is plotted in Fig. 3A, which shows that the CRSS of slip mode II increases with decreasing λ. When λ is less than ~15 nm, the CRSS of slip mode II is higher than that of slip mode I.

### CRSS of slip mode III

For dislocation slip mode III, its slip plane is parallel to the twin plane. Dislocation motion is blocked only by grain boundaries but not by TBs. As a result, the CRSS is dependent on grain size but independent of twin thickness. To obtain the CRSS of slip mode III, a dislocation pileup model in the presence of grain boundaries is used. A schematic illustration is shown in the inset of Fig. 3B. The CRSS of this slip mode can be expressed as follows (*10*, *11*)(7)where *K* is the Hall-Petch slope due to decreasing grain size *d*. In slip mode III, dislocations can slip within the twin plane and in the grains’ interior. The corresponding lattice frictional stress for shuffle-set 0° perfect dislocation BC in the twin plane and the interior of a grain is 17.2 and 10.3 GPa, respectively. Therefore, τ_{0} = 10.3 GPa, which is equal to the lattice frictional stress for dislocation slip in crystal interior. The Hall-Petch slope *K* is calculated as follows (*28*)(8)and based on the parameters of diamond (table S3), the CRSS is(9)The calculated CRSS increases with decreasing *d* (Fig. 3B), showing excellent agreement with experimental results for polycrystalline diamond samples (*6*, *7*).

### Origin of the unprecedented hardness in nt-diamond

On the basis of the CRSSs for the three dislocation slip modes, yield strength and hardness of nt-diamond can be evaluated using the Sachs model (*29*). We have constructed nt-diamond models with 6000 gains (with average grain sizes of 20 and 125 nm) with random orientation and various twin thicknesses (inset of Fig. 4A) to calculate yield stresses and hardness by evaluating the uniaxial stress for 90% deformed grains. Note that the grain sizes of nt-diamond used in this model are greater than those of nanocrystalline diamond with highest hardness reported (*7*); reverse Hall-Petch effect is not considered here (*30*). As shown in Fig. 4A, with increasing uniaxial stress, the population of plastically deformed grains increases owing to different grain orientations and CRSSs. When the population of plastically deformed grains reaches 90% (*31*), the corresponding uniaxial stress can be considered as yield stress, which is approximately one-third of the hardness (*32*). The calculated hardness of nt-diamond is plotted in Fig. 4B and fig. S4. The hardness for the model with average grain size of 20 nm agrees well with the experimental results for nt-diamond with the same grain size (shown in Fig. 4B) (*6*–*9*). This suggests that the model and parameters used are reliable. The results also indicate that if grain boundaries are all locked (that is, no grain boundary sliding) and only the three slip modes are considered, then a maximum hardness of 395 GPa is expected for nt-diamond with a grain size of 20 nm and twin thickness of 0.62 nm. Although the hardness of nt-diamond monotonically increases with decreasing twin thickness, the increasing slope is different at various twin thickness ranges. This is attributable to different contributions of the slip modes at different twin thickness ranges. At a grain size of 125 nm, the hardness of nt-diamond (shown in fig. S4) is slightly less than that of the 20-nm grain size model at the same twin thickness. This is attributed to the CRSS for slip mode III increasing with decreasing of grain size according to the Hall-Petch effect.

For the three slip modes, all of their CRSS values consist of two parts: lattice frictional stress τ_{0} and enhanced stress (or athermal stress τ_{u}) owing to effects of TBs and/or grain boundaries. The ratio of τ_{u}/τ_{0} increases with decreasing twin thickness λ for nt-diamond with a grain size of 20 nm. When λ decreases from 10 to 0.618 nm, τ_{u}/τ_{0} increases from 0.8 to 3.3. This shows that contributions of lattice frictional stress and athermal stress to hardness of nt-diamond have similar magnitudes. For metals, the behavior is totally different. In nt-Cu (*26*), for example, the strength is mainly contributed by athermal stress. Therefore, the origin of the unprecedented hardness of nt-diamond mainly comes from two factors: high lattice frictional stress, which is the nature of diamond bonding, and high athermal stress, which is due to the Hall-Petch effect. The Hall-Petch effect, which may amount to more than three times the lattice friction stress, plays an overwhelming role in enhancing the hardness of nt-diamond. Because of the covalent bonding in diamond, dislocation nucleation and motion in slip planes have low activation volume and high activation energy, as shown in fig. S3, which leads to activation energy insensitivity on external stress which gives high lattice frictional and athermal stresses appear in nt-diamond (*33*).

In summary, three dislocation slip modes in nt-diamond—slip transfer mode, confined layer slip mode, and paralleled to twin plane slip mode—are identified by considering orientation of dislocation lines and their kink with respect to twin planes. The interaction between dislocations and twin planes is studied by investigating the associated reaction heat, activation energies, and barrier strength using the MD method. CRSSs for the three dislocation slip modes have been obtained. On the basis of the Sachs model, hardness of nt-diamond is calculated. The results agree well with experimental data. The unprecedented hardness in nt-diamond is mainly contributed by two parts: high lattice frictional stress, which is due to the covalent nature of diamond, and high athermal stress, which is due to the Hall-Petch effect. This work provides a new insight into nt-diamond hardening mechanism, which will be helpful for studying superhard materials in experiments.

## MATERIALS AND METHODS

### Computational setup for MD simulation

To study dislocation interactions with twin planes, we built a series of atomic structure models consisting of coherent *Σ*3 TBs with various dislocation types (as shown in Fig. 1C) (*34*, *35*). A typical atomic structure model is shown in Fig. 1C. It contains 115,200 atoms, and the *x*, *y*, and *z* axes are along the , , and [111] directions, respectively, of the diamond matrix. MD simulations were performed using the LAMMPS (Large-scale Atomic/Molecular Massively Parallel Simulator) program (*36*), and C–C bonding interactions were described by the LCBOP (Long-range Bond-order Potential for Carbon) potential (*37*). Periodic boundary condition was only imposed along the *y* axis, and free surface was imposed in the *x* and *z* directions (Fig. 1C). A screw dislocation line was constructed in the twin plane. An enlarged dislocation structure is shown in Fig. 1D. As shown in Fig. 1D, the dislocation line is along the *y* axis. The constructed structures were relaxed via energy minimization.

### Computational method for reaction heat of dislocation reaction

Reaction heat of dislocation TB reaction in slip mode I was obtained by calculating the total energy difference before and after the dislocation reaction. To exclude the TB and the surface effect on the dislocation energy, the dislocation energy was calculated by setting the dislocation line as far as possible from the TB and free surface. An energy convergent test was performed by using a serial of supercells with different sizes. The calculated dislocation energy as a function of supercell size is plotted in fig. S1. When the supercell size is larger than 198 Å × 25 Å × 98 Å, the calculated dislocation energy will converge. Therefore, the supercell size 209 Å × 25 Å × 138 Å (contains 115,200 atoms) was used, which is larger enough to exclude effect of boundaries.

### Computational method for activation energy of dislocation reaction

The process of dislocation reaction with twin planes is a process of dislocation kink nucleation and migration in essence. Therefore, the activation energy of dislocation reaction with the twin plane is the activation energy for dislocation kink nucleation and migration. Schematics of the kink motion path for the dislocation interacting with the twin plane are plotted in fig. S2. For a given dislocation reaction with the twin plane, several stable configurations of dislocations with different kink pair widths were constructed. By using the nudged elastic band method, the activation energy for dislocation kink nucleation (*E*_{f}) and migration (*E*_{m}) can be obtained (*38*). The activation energy *Q* of the dislocation reaction with the twin plane was then obtained according to *Q* = 2*E*_{f} + *E*_{m}(*39*). The calculated results are summarized in Table 1 and table S1.

### Computational method for barrier strength of dislocation reaction

Barrier strength is defined as the threshold stress of dislocation reaction with TB when the activation energy is equal to zero. To obtain the barrier strength for the dislocation reaction with twin plane, dislocation reaction activation energy was calculated by adding a serial of shear stress to the MD simulation cell by using the method as described above (computational method for activation energy of dislocation reaction). Finally, the resolved shear stress–dependent activation energy of dislocation interaction with twin plane can be obtained. The obtained shear stress–dependent activation energy is fitted by using the following expression (*38*)(10)where *Q*_{0} is the activation energy at temperature of 0 K and stress of 0 GPa; τ_{TB} is the barrier strength; *p* and *q* are the energy barrier shape parameters. The fitted parameters in Eq. 10 for different dislocation reactions are plotted in fig. S3 and summarized in table S2. When the activation energy reaches zero, the corresponding shear stress is the barrier strength.

## SUPPLEMENTARY MATERIALS

Supplementary material for this article is available at http://advances.sciencemag.org/cgi/content/full/4/9/eaat8195/DC1

Fig. S1. Dislocation energy as function of supercell size.

Fig. S2. Schematic of process for dislocation reaction with twin plane.

Fig. S3. Calculated shear stress–dependent activation energy of dislocation reaction with twin plane.

Fig. S4. Hardness of nt-diamond at grain size of 125 nm.

Table S1. Calculated kink energy *E*_{f} and activation energy *W*_{m} of kink migration for dislocation reactions by slip transfer mode.

Table S2. Fitted parameters used in Eq. 10.

Table S3. The parameters used to calculate the critical resolved shear stress of slip transfer mode, confined layer slip mode, and paralleled to twin plane slip mode.

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is **not** for commercial advantage and provided the original work is properly cited.

## REFERENCES AND NOTES

**Acknowledgments:**We thank two anonymous reviewers for their feedback on this manuscript.

**Funding:**This work supported by the National Natural Science Foundation of China (grant nos. 51771165, 51372215, and 51332005) and NSF (grant no. EAR-1361276).

**Author contributions:**B.W. conceived the project. J.X., H.Y., and X.W. constructed the dislocation structure. J.X. performed all calculations. J.X., B.W., F.Y., and P.L. analyzed the calculated results. J.X., B.W., Y.W., X.Z., and Y.T. co-wrote the paper. All authors discussed the results and commented on the manuscript.

**Competing interests:**The authors declare that they have no competing interests.

**Data and materials availability:**All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.

- Copyright © 2018 The Authors, some rights reserved; exclusive licensee American Association for the Advancement of Science. No claim to original U.S. Government Works. Distributed under a Creative Commons Attribution NonCommercial License 4.0 (CC BY-NC).