Directional bonding explains high conductance values of atomic contacts in bcc metals

TL;DR

This paper demonstrates that directional bonding in bcc metals leads to highly coordinated stable atomic contacts with high conductance, explaining experimental observations that contradict traditional expectations.

Contribution

It introduces a classical molecular dynamics approach with directional bonding potentials to explain high conductance in bcc metal contacts, highlighting the importance of bonding directionality.

Findings

Highly coordinated structures form before rupture due to directional bonding.

Conductance calculations match experimental measurements for bcc metals.

Directional bonding significantly influences atomic contact stability and conductance.

Abstract

Atomic-sized junctions of iron, created by controlled rupture, present unusually high values of conductance compared to other metals. This result is counter-intuitive since, at the nanoscale, body-centered cubic metals are expected to exhibit lower coordination than face-centered cubic metals. In this work, classical molecular dynamics simulations of contact rupture, using an interatomic potential that accounts for directional bonding, yield highly-coordinated stable structures before rupture, unlike an isotropic bonding potential, which results in the expected stable single-atom contacts. Density functional theory electronic transport calculations show that conductance values of these highly coordinated and highly stable structures, can explain the experimentally measured values for conductance of body-centered cubic atomic contacts, thus revealing the important role of directional bonding in these metals.

Equations

$$F=\frac{{\Sigma }_{n,\sigma }{T}_{n,\sigma }(1-T_{n,\sigma })}{{\Sigma }_{n,\sigma }{T}_{n,\sigma }}$$

Full text

††preprint: APS/123-QED

Stretching a metallic nanowire results in a progressive reduction of its cross-section at the weakest point, until it finally breaks. From an atomistic viewpoint, when the minimum cross-section of the nanowire contains only a few atoms, and for very slow stretching, the minimum cross-section can in fact decrease by one atom at a time Landman et al. (1990); Agraït et al. (2003). It seems reasonable to assume that the ultimate stable contact that holds the metal together is a single atom. Indeed, measurements of conductance for atomically sharp contacts seem to point in this direction since the stable contact conductance before rupture for most metals is just above one quantum, except for some notable exceptions; such as, iron Vardimon et al. (2013); Calvo et al. (2009), tantalum, molybdenum and tungsten den Boer et al. (2007); Halbritter et al. (2003). In fact, it is still widely assumed that the chemical valence of the bridging atom in single-atom contacts primarily determines the pre-rupture conductance values Requist et al. (2016); Scheer et al. (1998).

Atomic-sized contacts are typically realized via mechanically controllable break junctions (MCBJ) Muller et al. (1992a, b) or scanning tunneling microscope break junctions (STM-BJ) Pascual et al. (1993); Agraït et al. (1993) and characterized through electron transport measurements. Classical molecular dynamics (CMD) simulations and first principles transport calculations have been key in providing interpretations of the experimental results Landman et al. (1990). In the past decade, through the combination of experiments and simulations, researchers have identified different atomic contact geometries that may form just before rupture Hasmy et al. (2005, 2001); García-Mochales et al. (2005); Sabater et al. (2013) and stated the important role played by the nearest neighbor atoms in the mechanical and electrical properties of these structures Sabater et al. (2018).

Computationally, most of the systems studied so far have been face-centered cubic (FCC) crystals, which exhibit a high probability to form single-atom contacts immediately before rupture. However, recent calculations Vardimon et al. (2016) show a clear mismatch between experiment and theory in the case of body-centered cubic (BCC) iron. The reason behind this discrepancy is not clear, since three factors could play a role: chemical valence, crystal lattice structure and/or magnetism.

To shed light on this question, we first perform CMD simulations of the iron rupture process, using two different interatomic potentials: one in which the bonding between atoms is treated as isotropic, the other in which the bonds also have covalent character. Density functional theory (DFT) calculations on CMD snapshots of the atomic configurations are then used to obtain the electronic transport properties (the conductance) and to compare with experimental data obtained from electron transport experiments in an STM-BJ at low temperature (4.2 K).

In the STM-BJ configuration, the electrode tip can collide with and be withdrawn from the surface over continuous cycles of rupture and formation of the atomic-sized contact. The electrode tip consists of iron wire of nominal diameter 0.25 mm (with a purity of 99.99%), and is connected in an electrical circuit as shown in Fig. 1 a), where a constant bias voltage $V=100$ mV is connected in series with the IV converter amplifier and, in turn, to a resistor of interest (in our case an atomic-sized contact). We express conductance that is the inverse of the resistance, in units of the quantum of conductance $G_{0}=2e^{2}/h$, where $e$ is the charge of the electron, $h$ is the Planck constant and the factor of 2 accounts for the spin degeneracy.

Typically, we record the conductance as a function of the relative displacement between the electrodes. The resulting curves are called rupture traces, as shown in Fig. 1 b). From each trace, we can build a histogram of conductance, which reveals the most frequent conductance values that the trace contains. Upon accumulating a significant number of individual histograms (one for every trace), we can construct a full histogram of conductance. For iron, it exhibits a clear peak at $\approx 2G_{0}$ (see, Fig. 1 c)), as has been previously observed Untiedt et al. (2004); Vardimon et al. (2016); Calvo et al. (2009).

In CMD simulations, the trajectory of each atom in the contact can be obtained by solving Newton’s second law by using a suitable interatomic potential to describe the interactions between the atoms. This constitutes the basic principle of classical molecular dynamics Rapaport (2004), which provides us with the means to model a very wide variety of materials down to the atomic level. The real problem here, however, is to select the most appropriate semi-empirical interatomic potential Daw and Baskes (1983), such that the interactions between the iron atoms are described with as much detail as possible. By far the most extensively used empirical potential to model metals is the embedded-atom method (EAM) potential Malerba et al. (2010); Sabater et al. (2018, 2013); Vardimon et al. (2016); Cuevas and Scheer (2010); Agraït et al. (2003), in which bonding is assumed to be isotropic. While this assumption is valid for FCC structures, BCC structures have a lower coordination about an individual atom and thus exhibit more directionality in their bonding, i.e. slightly covalent character Pettifor (1996). Accordingly, the modified embedded-atom method (MEAM) potential Baskes (1992) may provide a more realistic description of the bonding in BCC metals, since directionality is included. In this work we thus compare the MEAM and EAM methods in order to determine which will lead to better agreement with experimental results in electronic transport calculations on CMD-generated structures.

For both potentials, we use the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) Plimpton (1995); Plimpton et al. (2016). Additionally, to imitate the experimental conditions, all the simulations are realized using a Nosé-Hoover thermostat Nosé (1984); Hoover (1985) to maintain a constant temperature. The thermostat is applied at the recommended interval of 1000 simulation time steps Plimpton et al. (2016). We use a time step of 1 fs and the same initial input structure, consisting of $\approx 1500$ atoms, for comparison of the two potentials. Figure 2 shows a representative example of the rupture process of iron using the MEAM potential, with the atoms initially occupying positions in a perfect BCC lattice oriented along the (001) crystallographic direction (see Fig. 2 a)). The initial velocities of the atoms are randomized at the beginning of each rupture run and correspond to an average temperature of 4.2 K. The input structure is stretched at $\approx 1$ m/s until rupture. During every single rupture simulation, out of an ensemble of 100 independent runs performed with each potential, we compute the number of atoms in the minimum cross-section of the model contact by means of the Bratkovksy algorithm Bratkovsky et al. (1995). The minimum cross-section and simulation trajectory are both recorded every picosecond. For the purpose of comparison, traces are truncated 100 ps before the moment of rupture when constructing cross-section histograms.

Figure 3 compares the two normalized histograms of minimum cross-section data obtained with each of the two potentials. The pink-shaded histogram in Fig. 3 has been constructed by using the most recent MEAM interatomic potential, fitted to the melting point of Fe as well as its near-melting point elastic constants Etesami and Asadi (2018). This particular potential is suitable for simulations of Fe contact rupture because the (001), (110) and (111) exposed surface energies agrees very well with experiments Etesami and Asadi (2018). For comparison, see the blue-shaded minimum cross-section histogram in Fig. 3, obtained using an EAM potential whose surface energies also agree reasonably well with experiment and DFT calculations Malerba et al. (2010). Although both potentials show the formation of one and two-atom contacts (first two peaks in the histogram), the MEAM potential produces stable pre-rupture structures with a higher number of atoms in the minimum cross-section than the EAM potential. This is shown in the shaded tail of the histogram, where the probability of forming structures with more than two atoms right before rupture is lower with the EAM potential than with MEAM.

Note that in the case of FCC metals Hasmy et al. (2001, 2005); García-Mochales et al. (2005); Sabater et al. (2018), and even in the structures obtained with the EAM potential for Fe, the simulated cross-section narrows atom by atom during the rupture process. In contrast, the cross-section obtained using the MEAM potential breaks through cleavage across (110) oriented planes as shown in Fig. 2 e). These (110) planes are formed after a re-orientation of the contact during the applied tension.

To obtain the conductance of snapshots extracted from CMD simulations, such as in Fig. 2 c), we use the electronic transport code Alicante nanotransport (ANT.Gaussian) Palacios et al. (2001, 2002); Louis et al. (2003); Palacios et al. ; Frisch et al. . (For more details of the DFT calculations, see the Supplemental Material SM .) Conductance values near the peak of the experimental histogram in Fig. 1 c) are obtained for structures which correspond to those shown in Fig. 2 c) (see cases marked with an asterisk (*) in table SI of the Sup. Mat. SM ), that is, those predicted by the MEAM potential. Note that the EAM potential only reproduces the rupture process illustrated in Fig. 2 in 3 out of the 100 rupture simulations versus 20 out of a 100 in the case of the MEAM potential.

More revealing than the statistical study presented above, are the so-called Fano factors Cuevas and Scheer (2010). In experiments at low temperatures, low bias voltage and low frequency range (maximum of 400 kHz), the Fano factor $F$ provides a measure of noise suppression relative to the maximum Poissonian value of $2eI$ Cuevas and Scheer (2010). Therefore, the shot noise from ballistic transport of an electron through an atomic-sized contact, is given by $S_{I}=2eIF$, where $I$ is the bias current and $e$ the electron charge. On the other hand, in DFT quantum transport calculations, the spin-polarized conductance can be expressed as Cuevas and Scheer (2010): $G=\frac{G_{0}}{2}\Sigma_{n,\sigma}T_{n,\sigma}$, where $G_{0}$ is the usual spin-degenerate quantum of conductance and $T_{n,\sigma}$ are the individual spin-resolved eigenchannel transmissions Jacob and Palacios (2006). Since not only the geometry but also the number of atoms in the constriction of a contact determines the overall conductance through their valence orbitals, the individual spin-resolved transmission channels can convey information about the atomic structure of the contacts through the Fano factor (a measure of the number of partially open transmission channels in an atomic-sized contact):

[TABLE]

Figure 4 shows an example of how an eigenchannel analysis can be carried out to obtain the Fano factor from a conductance calculation on a model contact. A conductance calculation on the structure shown in Fig. 4 a) not only yields the overall spin-resolved transmissions in b) (which sum to $G=1.6574G_{0}$), but also the contributing spin-resolved eigenchannel transmissions, shown in Fig. 4 c). Therefore, in this case only 5 spin-resolved eigenchannels –3 spin-majority (purple) and 2 spin-minority (blue)– contribute significantly to the overall transmission. Based on this analysis, one can conclude that at least 2 atoms effectively contribute to the transmission in this contact. The noise in the transmission functions in Fig. 4 b) result from the disorder in CMD structures in general Calvo et al. (2008) and, in particular, from the spd hybridization of the spin-minority channels García et al. (2000).

Fano factors calculated from CMD snapshots of the two CMD potentials are presented in Fig. 5. The dark grey areas delineate forbidden values of $F$ for magnetic atomic-sized contacts, while the areas underneath the light grey line are the forbidden values of $F$ for non-magnetic materials. Recall that $F$ is a measure of the number of partially open transmission channels in a contact, and the more channels contribute to the overall conductance, the more atoms are likely involved.

For the MEAM potential (Fig. 5 b)), the calculated conductance values fall at or near the 4, 5 and 6 transmission channel lines. Experimentally, Fe has been shown to form last-contact structures with 6 transmission channels (see the experimental Fano diagram in Fig. 6 b) of ref. Vardimon et al. (2016)), which appears to indicate the formation of contacts with 3 or more atoms in them, according to the values of the Bratkovksy minimum cross-sections in table SI of the Supplemental Material SM . The MEAM potential, with more covalent character, is thus seen to outperform the EAM potential in this regard, whose $F$ values for the 17 conductance values collected in table SII SM are plotted in Fig. 5 a).

Moreover, based on the low density of states of Fe at the Fermi level García et al. (2000), in comparison to Ni or Co, one would expect Fe to have a first maximum conductance peak at a lower conductance value than either of the latter metals. The experimental Fano diagram in Fig. 6 b) of ref. Vardimon et al. (2016) exhibits a significant number of conductance values at this expected low value of $\approx 1.2-1.4G_{0}$, but in a histogram, they are subsumed by the broad peak at $\approx 2G_{0}$ (see fits to the histogram of conductance in the Supplemental Material SM ). Our interpretation, based on the simulation results presented above, is that slight differences in the structures responsible for this peak (see Fig. 2) could result also in deviations from the $2G_{0}$ value.

Therefore, we postulate that the discrepancy between the experiments and the simulations in the work of Vardimon et al. Vardimon et al. (2016) for the case of Fe is the lack of CMD structures with minimum cross-sections above $\approx 3$ atoms at rupture when using an EAM potential, in a combination with a tight-binding model to calculate conductance. Using a more detailed model of the bonding between Fe atoms, such as the MEAM potential used in this work, should markedly improve the comparison between the experimental histogram of conductance and the one calculated in Vardimon et al. (2016) based on the good agreement between the theoretical Fano diagram in Fig. 5 b) and their experimental Fano diagram in Fig. 6 b) of Vardimon et al. (2016).

It is also important to note that other BCC materials such as Ta, Mo and W, in similarity with iron, show experimental histograms of conductance with pronounced peaks at around $\approx 2G_{0}$ den Boer et al. (2007); Halbritter et al. (2003), and that the formation of similar structures could be explored in future work.

In summary, the pronounced peaks at $\approx 2G_{0}$ that appear in experimental histograms of conductance of body-centered cubic metal atomic-contacts could not be explained by considering single-atom contact structures Etesami and Asadi (2018). Here, we show that DFT electronic transport calculations of structures with several atoms in the minimum cross-section give conductance values in agreement with experiments. Such structures arise when an energetically favorable reorientation occurs, from (001) to stable (011) layers perpendicular to the direction of stretching. We obtain this stable reorientation when using an interatomic potential (the MEAM potential) that includes directional bonding. Our findings contradict the presently-held assumption that the most likely stable pre-rupture contacts in BCC metals are made up of a single atom. The MEAM potential thus provides a realistic mechanism of atomic rupture for Fe in which covalent bonding plays a key role. Body-centered cubic metals such as Fe may therefore represent good candidates for producing molecular junctions in which the electrode-molecule interface is an atomically flat surface.

I Acknowledgments

This work was supported by the Generalitat Valenciana through PROMETEO2017/139 and GENT (CDEIGENT2018/028), the Spanish government through grants MAT2016-78625-C2-1-P, FIS2016-80434-P and the María de Maeztu Programme for Units of Excellence in R&D (MDM-2014-0377), by Comunidad Autónoma de Madrid through Grant S2018/NMT-4321 (NanomagCOST-CM), by the Fundación Ramón Areces, and by the European Union Graphene Flagship under Grant No. 604391. JJP acknowledges the computer resources and assistance provided by the Centro de Computación Científica of the Universidad Autónoma de Madrid and the RES.