Method for determining optimal supercell representation of interfaces

The geometry and structure of an interface ultimately determines the behavior of devices at the nanoscale. We present a generic method to determine the possible lattice matches between two arbitrary surfaces and to calculate the strain of the corresponding matched interface. We apply this method to explore two relevant classes of interfaces for which accurate structural measurements of the interface are available: (i) the interface between pentacene crystals and the (111) surface of gold, and (ii) the interface between the semiconductor indium-arsenide and aluminum. For both systems, we demonstrate that the presented method predicts interface geometries in good agreement with those measured experimentally, which present nontrivial matching characteristics and would be difficult to guess without relying on automated structure-searching methods.

The geometry and structure of an interface ultimately determines the behavior of devices at the nanoscale. We present a generic method to determine the possible lattice matches between two arbitrary surfaces and to calculate the strain of the corresponding matched interface. We apply this method to explore two relevant classes of interfaces for which accurate structural measurements of the interface are available: (i) the interface between pentacene crystals and the (111) surface of gold, and (ii) the interface between the semiconductor indium-arsenide and aluminum. For both systems, we demonstrate that the presented method predicts interface geometries in good agreement with those measured experimentally, which present nontrivial matching characteristics and would be difficult to guess without relying on automated structure-searching methods.

I. INTRODUCTION
As electronic devices shrink in size to reach nanoscale dimensions, interfaces between different materials become increasingly important in defining the device characteristics 1 . In many cases, it has been shown that the effect of the interface even dominates the device properties 2 , leading to the concept that "the interface is the device" 3 . In order to optimize the performance of a device it is therefore important to understand the properties of its interfaces.
First principles modeling based on atomistic methods such as Density Functional Theory (DFT) 4 have become an important tool for simulating the properties of interfaces 5 . To be truly predictive, atomistic methods require an accurate model for the atomic-scale geometry of the interface. As these simulations typically use periodic boundary conditions in the direction parallel to the interface, a common supercell for the surfaces of the two crystals forming the interface must be determined. However, typically the two crystals are not commensurate and finding a common supercell requires straining one of or both the surfaces. To accommodate the resulting strain, the two surfaces can also be rotated with respect to each other. However, for rotation angles preserving a high symmetry in the supercell, this has often the side-effect of increasing considerably its dimensions. Finding a supercell with low built-in strain and without an excessive number of atoms is therefore highly nontrivial.
In this paper we present an algorithm which allows for an efficient and systematic search for common supercells between two crystalline surfaces. Given the optimized geometries of two surfaces forming the interface, the algorithm returns a list of all possible interface supercells by varying the interface strain and the rotation between the two surfaces. A related, but more simplistic method have been proposed in Ref. 6 . Compared to Ref. 6 our method automatically tests all possible rotations of the two lattices and has been implemented into a graphical user interface, the Virtual NanoLab 7 .
In the paper we show that this is not only a practical procedure for generating low strain supercells for atomicscale simulations, but is a predictive tool for determin-ing interface geometries in accordance with experimental data. As a first example, we consider the interface between a pentacene crystal (PC) and the Au(111) surface, which has been widely studied both theoretically [8][9][10][11][12][13] and experimentally [14][15][16][17][18][19][20][21][22][23][24] . We show that the predicted geometries of a pentacene monolayer on Au(111) recover those observed experimentally. Using DFT, we calculate the ground state structure and energetics of these interfaces and find that they are thermodynamically more stable than those previously used in the literature.
As a second example, we consider the interface between Al and InAs. This interface is relevant for studies on semiconductor nanowires (NWs) in which superconducting properties are introduced by proximity effect with a superconductor 25-28 and its structure has been recently resolved using high-resolution transmission electron microscopy (HR-TEM) 29 .
The organisation of the paper is the following. In section II we introduce the algorithm for matching the two crystal orientations with minimal strain. In section III A we first apply the method to determine the geometry of a pentacene overlayer on Au(111). In section III B we determine the structure of Al on InAs. Finally in section IV we conclude.

II. METHODS
Algorithm for surface matching: A scheme of our algorithm for surface matching is shown in Fig. 1. In order to obtain a low-strain interface structure, we systematically search through all possible 2D unit cells of the two surfaces forming the interface. Given two arbitrary surfaces A and B with primitive vectors ( a 1 , a 2 ) and ( b 1 , b 2 ), we generate the Bravais lattices of the possible surface supercells A * and B * :  For each pair of supercells, we next determine a rotation matrix R which rotates B * and aligns u 1 with v 1 : , respectively. Finally, we match the two supercells by defining a strain tensor , which is applied to B * in order to match its Bravais lattice onto that of A * . The resulting equation to match A * and B * reads: with the individual components of the strain tensor ε being: A similar procedure can also be be applied to strain A * and match it to B * , or to strain equally both surfaces.
Straining one or both sides of the interface introduces and additional elastic contribution to the interface energetics. This contribution and its influence on the surface geometry varies considerably depending on the strength of the interaction between the overlayer and the substrate. For a substrate and a strained overlayer with cubic symmetry, we can write down the total energy per unit of area as 30 : where E int is the interface energy between the substrate and the overlayer, E surf is the energy of the free surface of the overlayer, ε xx , ε yy , ε xy have been defined in Eq. 5-7, C 11 , C 12 , C 44 are the elastic constants of the overlayer material and t is its thickness. The equation can be further simplified as Neglecting interactions between the interface and the overlayer free surface and strain effects, E surf will be independent of the interface geometry. For metals and weakly interacting interfaces, we also expect the interface energy E int to be rather similar for different geometries, so that the contribution of the elastic energy will be dominant and will determine the stability trend of the different geometries. On the other hand, for interfaces between semiconductors there may be a varying number of bonds at the interface, depending on the precise overlay structure. Since binding energies for covalent bonds are typically in the range 1-2 eV, the overlayer may need to have a thickness above ∼2 nm before the elastic energy will dominate.
We have implemented the algorithm into the Virtual NanoLab 7 and the calculated matches are presented graphically as illustrated in Fig. 3, The plot shows the number of atoms in the supercell as function of the average strainε. For the average strain we for simplicity use,ε = (|ε 11 | + |ε 22 | + |ε 12 |)/3. The measure in Eq. (10) gives slightly different orderings, but we found that orderings are basically similar and therefore selected the most simple option. The algorithm has been implemented to only test relevant vectors and with default values of N max , M max = 6 a scan typically takes 1s and the simulation time scale as N 2 max . In the next section we will apply the method to determine the structure of a Pentacene monolayer on Au(111) and the interface geometry of Al on top of InAs.

A. PC/Au(111)
Computational details: The DFT calculations for PC/Au(111) have been performed using the Atomistix ToolKit 31 . The Kohn-Sham orbitals have been ex-panded in a linear combination of pseudo-atomic orbitals (PAOs) 32 . The electronic exchange-correlation (xc) energy has been described by using the generalized gradient approximation (GGA) and the Perdew-Burke 91 (PW91) xc-functional 33 . We use this functional to compare with previous calculations 9 , however, note that the GGA-PW91 xc-functional does not include van-der-Waals forces and it will therefore underestimate adsorption energies. We have used a slab geometry with periodic boundary conditions parallel to the surface and mixed (Dirichlet+Neumann) boundary conditions in the direction normal to the surface, the latter allows for describing slabs with different workfunctions on the upper and lower surface. The Brillouin zone has been sampled using an 8 × 3 × 1 Monkhorst-Pack 34 grid and a Fermi-Dirac occupation scheme with a broadening of k B T = 25 meV. Structural relaxations have been performed using a convergence threshold for the forces of 0.01 eV/Å. During both the structural optimization and the evaluation of binding energies, the basis set superimposition error (BSSE) has been corrected using the counterpoise (CP) correction scheme 35 . For the Au(111) surface, only the two uppermost layers were allowed to relax during the structural optimizations, while the atoms in the lowermost layers were kept frozen at their bulk position.
For carbon we have used 21 orbitals per atom with s, p and d character and ranges up to 3.9Å, while for hydrogen we have used 5 orbitals per atom with s and p character and ranges up to 4.2Å. This basis set has been optimized to reproduce hydrogen and carbon dimer total energies 36 . Using this basis set, we obtain an adiabatic ionization energy for the individual pentacene molecule (P1) E I = 6.34 eV , in good agreement with the experimentally reported value E exp I = 6.59 eV 37 . For gold, we have used an s, p, d basis set of ranges 2.7 − 3.6Å, with a total of 9 orbitals per atom. The calculated lattice constant for bulk Au using this basis is a Au = 4.17Å. Using a layer of gold ghost orbitals to get a better description of the isolated Au(111) surface 38 , we also obtain that the surface work function is W Au(111) = 5.19 eV. Both values are in good agreement with those obtained using similar computational parameters and a plane wave basis set 9 .
Results: To construct the interface between the PC (see Fig.2) and Au(111), we have considered the PC unit cell according to the crystallographic parameters of Ref. We have then aligned the 010 direction of the PC unit cell along the normal to the Au(111) surface, and generated all possible interface supercells with N max , M max ≤ 12 by straining the PC lattice in the plane perpendicular to the 010 direction. In this case, we have used the experimental lattice parameter of the Au(111) surface. Fig. 3(a) shows a graph with the resulting possible supercells, sorted according to the mean absolute strain¯ and the number of atoms in the supercell. The lattices TABLE I. Strain in the 010 -oriented PC crystal to match Au(111). The first and second columns label each geometry by a roman number and list the supercell in the basis of the Au(111) bravais lattice. The third column list the number of PC(010) surface cells in the structure. ε11, ε22, ε12 are the components of the strain tensor applied to the PC(010) surface cell in order to match the gold supercell.ε = (|ε11| + |ε22| + |ε12|)/3, is the average strain.

Structure
Au ( with lowest strain are listed in Table I. It can be seen that supercell I, which has been used in earlier reports to model the PC/Au(111) surface 9 , possesses a rather large internal strain. On the other hand, our method reveals the existence of other non-trivial supercell arrangements which are associated with a much less strained PC lattice. In particular, among all the interfaces formed by a single PC(010) surface unit cell, the value of¯ is 67% and 77% lower for supercell II and III, respectively.
To analyze the relationship between strain and adsorption properties in PC/Au(111), we have compared the optimized geometries of supercells I, II, and III (see Table II). For each optimized geometry we calculate the lattice vectors a, b of the PC/Au(111) supercell, and the geometrical parameters of the PC crystal: the adsorption height z and the polar and azimuthal adsorption angles θ and φ, see Fig. 4(b,c). Since supercell II and III have very similar properties, in the following we will only compare supercell I and II, see Fig. 4 (a,b). The structural parameters obtained for supercell I are very similar to those obtained using a plane wave basis set 9 . However, supercell II and III provides an overall better agreement with the available experimental data. In particular, the supercell lattice vectors a, b, and the azimuthal angle φ are closer to those measured experimentally. In addition to the geometrical properties, we find that the calculated work function for supercell II is also in closer agreement with that measured experimentally, compared to that calculated for supercell I. Finally, the binding energy E b calculated for supercell II is also larger compared to that calculated for supercell I. This indicates that the supercell II and III, in addition to providing structural parameters in better agreement with those measured experimentally, lead to a structure which is thermodynamically more favorable. We note that the discrepancy relative to the experimental value, is due to the neglect of van der Waals forces.

B. Al/InAs
Following recent experimental work on InAs NWs in which superconducting properties have been induced by the proximity effect with Al 25-28 , we have considered two surfaces of InAs: the (1100) surface of the wurtzite phase (hereafter, (1100) WZ ), and the (111)B surface of the zinc-blend phase (hereafter, (111)B ZB ). NWs with both these surfaces orientations have been grown experimentally, and it has been demonstrated that the precise orientation of the epitaxial Al overlayer depends on the exposed InAs surface 29 .
For each of the two InAs surfaces, we have performed a scan over all Al(mkl) surfaces with m, k, l ≤ 3. Subse- For both the InAs surfaces considered, we have found that the Al surface which is predicted to have the lowest strain in the Al overlayer matches with that identified experimentally for thick (t > 30 nm) Al overlayers, see Table III. In the case of InAs(1100) WZ , this corresponds to Al(112), whereas in the case of InAs(111)B ZB , this corresponds to Al(111).
Another check on the accuracy of the method can be done by comparing the structures of the predicted InAs/Al interfaces with those measured experimentally. Fig. 5 shows the predicted InAs/Al interfaces, superimposed to the measured HR-TEM images for each interface. It can be seen that, for both interfaces, the agreement between the structural model and the HR-TEM pattern is excellent. On the InAs side of the interfaces, regions of dark and bright contrast can be associated with In and As atoms, respectively, whereas on the Al side of the interfaces, the regions with bright contrast surrounded by a darker halo can be associated with Al atoms.

IV. CONCLUSIONS
In conclusion, we have presented a systematic and efficient method for determining a supercell geometry of the interface between two crystals. The method has been implemented into the Virtual NanoLab software. The method was applied to two metal-semiconductor systems, Au-Pentacene and InAs-Al interfaces. In both cases the method suggests interface geometries in good agreement with experimental data. For Au-Pentacene we illustrated that previous studies 9 , which does not use a systematic approach for finding a supercell geometry of the interface have lower binding energies and are not in accordance with experimental data.