Quick-start guide for first-principles modelling of semiconductor interfaces

Interfaces between dissimilar materials control the transport of energy in a range of technologies including solar cells (electron transport), batteries (ion transport), and thermoelectrics (heat transport). Advances in computer power and algorithms means that first-principles models of interfacial processes in realistic systems are now possible using accurate approaches such as density functional theory. In this `quick-start guide', we discuss the best practice in how to construct atomic models between two materials and analysis techniques appropriate to probe changes in local bonding and electronic band offsets. A number of examples are given related to perovskite solar cells.

• mechanical strain, altering the local electronic structure; • chemical bonding environment, resulting in new electronic states active for trapping electrons and/or ions; • accumulation of intrinsic point defects and impurities, with enhanced mobility of defect centres.
The development of microscopic models for materials interfaces has a rich history. For example, in 1935 Gurney studied the effect of alkaline earth metals on the workfunction of tungsten using quantum mechanics 2 . Knowledge of valence and conduction band offsets contributed to the development of rectifiers and transistors. 3 In the late 1970s, the first descriptions of interfaces based on self-consistent density functional theory (DFT) were reported 4,5 . Procedures for the calculation of electronic band offsets were suggested later, 6,7 and were widely applied to compound semiconductors in the 1990s. 8,9 Beyond semiconductor-semiconductor interfaces, a wide range of solid interfaces have been studied including those involving metals, semiconductors and insulators, partly motivated by their importance in transistor technology. 5,[10][11][12][13][14] The latest generation of materials employed in energy storage and conversion -including perovskites, MXenes, spinels, olivines -are multicomponent and adopt crystal structures more complex a) Electronic mail: a.walsh@imperial.ac.uk than face-centred-cubic (fcc) semiconductors. 15 Investigation of realistic interfaces in these technologies (e.g. thermoelectrics, batteries, power electronics) pose scientific challenges. [16][17][18][19] Here, we discuss how to approach modelling and analysing interfacial processes using firstprinciples computational approaches.

I. ISOSTRUCTURAL INTERFACE MODELS
A special case for building interface models is when both materials share the same structure type. An atomic model can be readily constructed along different crystal orientations by a simple transformation of the unit cell. For example, a cubic unit cell is reoriented along < 111 > arXiv:1808.00359v1 [cond-mat.mtrl-sci] 1 Aug 2018 using a rotation matrix of: An example interface for fcc semiconductors is illustrated in Figure 1. A particular orientation may be favoured thermodynamically, following the most stable crystal surfaces, or may be determined by the fabrication procedure of a device. Relating to thin-film solar cells, both CH 3 NH 3 PbI 3 and CH 3 NH 3 PbBr 3 adopt a perovskite-type structure with a low lattice constant mismatch ( ∆a a < 3%), which can be used to construct small (< 100 atom) supercell models of the (001) interface that are practical for DFT calculations. 20 For ionic materials, care is needed as a particular crystal orientation may carry an electric dipole moment. 21 The classification of polar surfaces by Tasker 22 is summarised in Figure 2. For example, the (100) surface of rocksalt is non-polar (Type 1), the (110) surface is quadrupolar (Type 2), and the (111) surface is polar (Type 3).
There may be motivation to study a specific interface orientation (e.g. from microscopy), but in general a nonpolar orientation should be chosen to avoid complications in practical calculations. Type 3 surfaces generate a macroscopic electric field that is proportional to the size of the supercell expansion. In real materials, such an electric field is removed by structural or chemical reconstructions, but in calculations it can cause poor convergence in properties with respect to slab thickness, and unphysical behaviour such as charge 'sloshing' between the two sides of the interface.

II. HETEROSTRUCTURAL INTERFACE MODELS
For interfaces involving two materials with different crystal structures, challenges arise owing to the choice in the relative orientations of the two crystals. Often, experimental information on the preferred morphology is not available, especially for newer technologies.
Zur and McGill 23 proposed an algorithm lattice match to identify suitable periodic reconstructions of an interface between two dissimilar materials. A search is performed of lattice directions and expansions, and they are reduced to their simplest form. The example of Si on Al 2 O 3 is shown in Fig. 3. The approach was recently implemented by Butler et al to predict lattice-matched electrical contacts suitable for perovskite solar cells, which identified a favourable match between CH 3 NH 3 PbI 3 and Cu 2 O. 24 Other procedures for lattice matching have also been reported. 25,26 While procedures such as lattice match can identify low-strain orientations, the underlying chemical bonding (atomic sites and bonding) is neglected. Coincident Site Lattice Theory addresses the atomistic detail. 27 At certain misorientation angles, there may be more coincident lattice sites (i.e. a common sublattice) in addition to the origin, and such configurations can be stabilized over others with a similar misorientation angle. One might expect that an interface model with a high planar density of coincidence sites to be most stable; however, this strongly depends on how the atoms are bonded at the interface. For example, in an ionic material cations prefer to be in contact with anions. 28 In a recent study of interface formation between perovskite-structured CsPbBr 3 and rocksalt-structured PbS -a system with potential applications in light emitting diodes -a number of high-coordination configurations were identified and their stability assessed. 29 The interface configuration where Pb of CsPbBr 3 is located on S of PbS, and Br of CsPbBr 3 is located on Pb of PbS, was found to be most feasible due to favourable electrostatic interactions that lower the interface energy.

III. INTERFACE ANALYSIS
Once a reliable atomic model of an interface between materials X and Y has been constructed, and calculations have been successfully performed, properties can then be analysed. While some application areas require specific analysis -such as the potential energy barriers for ion diffusion across an interface in solid-state batteries 30 -the points discussed here are intended to be general.

A. Interface energy
Concerning thermodynamic stability, similar to the calculation of a defect formation energy, 31 the interface energy E f (X/Y ) can be calculated from the general expression: where E tot (X/Y ) is the total energy of interface model, A is the area of the interface in the supercell model, and µ i is the chemical potential of atomic species i. Since materials X and/or Y are strained biaxially in the interface calculations, the reference chemical potential is usually obtained from the strained bulk materials, and it is also important to consider relaxation in the direction normal to the strain plane (i.e. the Poisson effect). The chemical potential can also be used to account for different growth conditions, including contributions from off-stoichiometry and impurities.
Owing to periodic boundary conditions, two interfaces are usually contained within a supercell. Both X/Y and Y /X interfaces must be identical to obtain a well-defined interface energy for a particular interface; this corresponds to a supercell model with inversion symmetry. If the two interfaces are not identical, one can employ a vacuum slab geometry to obtain the interface energy for each interface in turn, while properly accounting for the associated surface energy.
If the models of X and Y are stoichiometric, then the chemical potential can be removed and the interface energy expression is simplified to: where n X is number of stoichiometric units of X and E tot (X) is the total energy of bulk X. Values of interface energy typically range from 0-5 Jm −2 (0.3 eV/Å 2 ). [32][33][34] Otherwise, the interface stability can be inferred from the work associated with placing the surfaces of X and Y in contact: where W sep (X/Y ) is the work of separation, E slab (X) and E slab (X) are the total energy of isolated X and Y surface slab models. Therefore, negative W sep (X/Y ) represents a preference for X and Y to form an interface. W sep is readily obtained from calculations, but in experiment cleavage of an interface also results in elastic, reconstructive, and diffusive processes that act to lower the energy of the cleaved surfaces. The measured quantity is called the work of adhesion W ad and differs from W sep by varying degrees; comparison between the two quantities should be made with caution. W sep is relevant for mechanical properties of interfaces; however, for dynamical processes such as surface wetting differences between W sep and W ad can be important. 35

B. Structure change
Interface simulations can be used to probe changes in bond lengths and angles across an interface, simply by inspection of the atomic coordinates. Due to the abrupt change in coordination environments across an interface, the differences before and after geometry optimisation can be large.
In addition to the equilibrium structure, the impact of tensile/compressive strain can be also investigated. If material Y is grown pseudomorphically on material X, then Y experiences biaxial strain, depending on the difference in lattice constant. The physical properties (e.g. band gap, point defect formation, and ion diffusivity) will be affected accordingly. Crystal strain could be beneficial, e.g. as employed in Si electronics where the channel Si is intentionally strained to increase the carrier mobility. 36 Strain at an interface can be reduced by elemental intermixing, 37 e.g. by having a X 1−y Y y compositional gradient rather than an abrupt junction; however, this can be challenging to model using atomistic simulation techniques due to the large model systems required.

C. Electrostatic potential
Beyond the local structure, the behaviour of mobile charge (ions and electrons) is influenced by the changes in electrostatic potential towards an interface. An example for the PbS | CsPbBr 3 interface is shown in Fig. 4. The 3D electrostatic potential fluctuates rapidly because of the atomic potentials. To obtain plateaus useful for alignments, the potential can be averaged over the supercell (length ∆) along the interface normal direction (z) as follows: 6 Plateaus are obtained in bulk-like regions if the dipoles at each interface cancel out. 9 Tools such as Macrodensity 38 allow for a range of potential averaging methods to be tested.
The electronic band edge positions of the bulk materials that form interface can be obtained using an appropriate alignment procedure; a selection of different approaches can be found in Refs. 7-9, 39-44. There are similarities to calculating the workfunction (or ionisation potential) of an isolated surface. An example is shown in Figure 5 for the CH 3 NH 3 PbI 3 (X = Cl, Br, I) series. 20 The valence and conduction band offsets were calculated using an alignment procedure that employs deep-lying core states, with results equivalent to the raw electrostatic (Hartree) potential. From this analysis, the band offset class can be assigned to distinguish between the cases of charge confinement ('straddling', type I), charge separation ('staggered', type II), or a conductive junction ('broken gap', type III). We note that the natural band offset procedure described above is not suitable for observing local band gap narrowing or widening at an interface, which can be analysed by the local density of states (DOS).

D. Electronic structure
In standard calculations, the electronic DOS is summed over the entire simulation cell; however, it can be projected in real space along the interface axis to probe changes in chemical bonding towards the junction.
Analysis of the charge density difference (∆ρ = ρ (X|Y ) − ρ X − ρ Y ) can be used to show how electrons are redistributed upon interface formation. If alternating charge accumulation and depletion regions are shown along the contact area, this implies chemical bonding that enhances interface stability. On the other hand, if a charge accumulation (depletion) region is shown alone, this could imply cation-cation (anion-anion) repulsion across a polar interface.
The magnitude of the electronic band gap can be increased or decreased in interface calculations because of strain, as mentioned above. Metal-induced gap states (MIGS) can also be introduced in the semiconductor side. 40 The MIGS result in charge redistribution and modified offsets at the interface. 46 DFT simulations with spatially resolved DOS can be used to characterise MIGS at metal-semiconductor interfaces. 47 At metal-semiconductor interfaces, an additional consideration is the image-potential. 48,49 In 1983, Stoneham proposed that the image potential could be the dominant contribution to surface adhesion and metal support effects in catalysts, where metal ceramic interfaces are present. 50 Later ab initio calculations of metal-ceramic interfaces confirmed that the image-potential is a leading term in determining interface stability in such systems. 51 Interface properties can also be investigated using the non-equilibrium Green's function (NEGF) formalism, which is often used for electron transport. A non-periodic interface is achieved by introducing semi-infinite electron reservoirs on each side of the simulation cell. Another advantage of this approach is that current-voltage characteristics can be compared directly with experiment. 45

IV. CONVERGENCE OF CALCULATIONS
There are standard convergence criteria for electronic structure calculations of bulk solids, which can include changes in total energy, ionic forces, and structure parameters with respect to the k -point sampling and the quality of basis set. The main two additional factors for interface calculations are slab orientation and thickness, which need to be converged with respect to the property of interest.
The description of electron transfer and space charging effects may be particularly sensitive to the supercell size, which is typically much smaller than the physical screening length for most materials.
There is also the choice of electron exchange and correlation functional in DFT. While one functional may provide a reliable description of bulk properties of one material, at the interface the chemical bonding is perturbed and the electronic wavefunction can be more localised. An accurate description of interfacial bonding may require hybrid functionals or beyond, depending on the chemistry of the system being studied. 52 When calculating separate representations of the bulk, surface and interface systems, it is advisable to check that all simulation cells are in the same crystallographic orientation. For reasons of numerical accuracy, ensuring that all cells have the same internal parameters (e.g. grid density and k -point sampling) leads to smoother convergence of the terms in Equations 2 -4.
For example, if simulating an interface along the < 111 > direction, the transformation in Equation 1 should be applied to the bulk as well as to the interface model. Additionally, the vacuum spacing in the surface model should be chosen to be an integer supercell expansion of the bulk to enable to clean comparison between the surface and interface models.

V. CONCLUSION
Procedures for simulating bulk crystalline materials are well-established. The combination of two or more materials to form an interface raises additional challenges for materials modelling, but can yield valuable insights that are difficult to probe experimentally. We have shown how consideration of the surface terminations and latticematching conditions are useful tools to generate interface models. Furthermore, we discussed and provided introductory references on how to approach analysis of the physical properties including thermodynamic stability, bonding, and electronic structure. There are many application areas for these techniques across emerging energy technologies.