Electronic structure and optical properties of Sr$_2$IrO$_4$ under epitaxial strain

We study the modification of the electronic structure in the strong spin-orbit coupled Sr$_2$IrO$_4$ by epitaxial strain using density functional methods. Structural optimization shows that strain changes the internal structural parameters such as the Ir-O-Ir bond angle, which has an important effect on the band structure. An interesting prediction is the $\Gamma - $X crossover of the valence band maximum with strain, while the conduction minimum at M remains unchanged. This in turn suggests strong strain dependence of the transport properties for the hole doped system, but not when the system is electron-doped. Taking the measured value of the $\Gamma-X$ separation for the unstrained case, we predict the $\Gamma - $X crossover of the valence band maximum to occur for the tensile epitaxial strain $e_{xx} \approx 3\%$. A minimal tight-binding model within the $J_{\rm eff} = 1/2$ subspace is developed to describe the main features of the band structure. The optical absorption spectra under epitaxial strain are computed using density-functional theory, which explains the observed anisotropy in the optical spectra with the polarization of the incident light. We show that the optical transitions between the Ir (d) states, which are dipole forbidden, can be explained in terms of the admixture of Ir (p) orbitals with the Ir (d) bands.


I. INTRODUCTION
The 5d oxides such as Sr 2 IrO 4 (SIO) are of considerable current interest due to the presence of a large spinorbit coupling (SOC) which leads to many novel features such as the spin-orbit assisted J eff = 1/2 Mott insulator 1 and spin-orbital entangled electron states. It has been suggested that the spin-orbital entanglement induced by the strong SOC in these structures could make these materials hosts for several unconventional features such as the Kitaev model 2,3 , quantum spin Hall effect at room temperature, 4 or unconventional superconductivity 5,6 .
Strain is an important parameter for probing the nature of the electron states, and it can be induced by pressure or by epitaxial growth on lattice-matched substrates. Tuning of the band gap is an important aspect for functional manipulation for potential device applications.
There have been several studies of SIO under epitaxial strain condition, both from theory and experiments. Samples of SIO have been grown epitaxially on a number of substrates such as SrTiO 3 , LaAlO 3 , GdScO 3 , etc. [7][8][9] . Resistivity and optical absorption measurements on these structures have shown that the Mott-Hubbard gap is preserved under epitaxial strain, but its magnitude can be tuned by varying the strain. The changes in the electronic structure show up in the optical properties as red or blue shift of the optical absorption under strain condition 8,10 . Several theoretical works have also addressed the electronic structure of SIO under strain [11][12][13] . The optical properties were studied by Zhang et al. 11 and Kim et al. 13 . These calculations were limited to the low energy range (0 -2 eV) and, furthermore, the polarization dependence of the optical absorption has not been studied theoretically, even though experiments show a strong anisotropy for absorption with E plane vs. E ẑ . Thus, a full understanding of the electronic and optical properties is still missing.
In this paper, we study the effect of the epitaxial strain on the electronic structure and optical absorption of SIO from density-functional theory. We find an interesting Γ − X crossover of the valence band top under strain, which we explain from a tight-binding model, and also find the anisotropy in the absorption spectrum for light polarized along the plane vs. normal to the plane, in agreement with the experiments. The dipole-forbidden d-d transitions are explained in terms of admixture of the Ir p orbitals with the d bands.
The organization of the paper is as follows. In Section II, we describe the crystal structure, the method of calculation, and the results of the density-functional structural optimization of the crystal structure under strain. Density-functional results for the band structure and magnetic moments under strain are discussed in Sec-tion III. Section IV discusses the optical absorption spectrum under strain, and the results are summarized in Section V.

II. CRYSTAL STRUCTURE AND METHOD OF CALCULATION
Bulk Crystal Structure -The crystal structure of SIO with space group (142) I4 1 /acd consists of eight formula units in the √ 2a 0 × c unit cell 14 as shown in Fig. 1, with the lattice constants a = b = 5.497Å and c = 25.798Å. (An equivalent unit cell with the body-centered tetragonal lattice and four formula units in the basis may also be used.) Note that a = √ 2a 0 . The structure shows a staggered rotation of the IrO 6 octahedra about the caxis by the angle θ = 11.5 • . The unit cell has four IrO 6 planes stacked along the c-axis, each plane consisting of two octahedra, with staggered rotation angles and antiferromagnetic Ir moments. The electronic structure is essentially controlled by the Ir atoms with the 5d 5 configuration placed in the crystal field of the oxygen octahedra, with the Sr atoms playing a passive role of donating electrons to the system, described by the nominal chemical formula Sr 2− 2 Ir 4+ O 2− 4 . The structure can be thought of as IrO 2 layers separated from each other by two intervening SrO planes, making the structure quasitwo dimensional (2D). This means that the basic electronic structure can be modeled by a single plane of the Ir atoms in a minimal model, which we discuss later. Density-Functional Methods -Density functional theory (DFT) with linearized full potential muffin-tin orbital method (FP-LMTO) [15][16][17] was used to solve the Kohn-Sham equations within the local spin-density approximation for the exchange and correlational functional 18,19 . The LMTO basis set consisted of spdf orbitals for Ir and Sr and spd orbitals for O, which were augmented outside the muffin-tin spheres by two Hankel functions of two different decay lengths. The semicore states Ir 5p and Sr 3p were treated as valence electrons. The band calculations were carried out for the optimized structures within the local spin density approximation including both the SOC and the Hubbard U terms (LSDA + SO + U). Following earlier authors, we used U = 2.7 eV. We also employed Vienna ab initio simulation package (VASP) in projector augmented wave (PAW) formalism including Hubbard U = 2 eV and spin-orbit interaction 20 for computing magnetic moments.
Structural optimization -Atomic positions were relaxed in all our calculations within the local spin-density approximation and the force convergence was obtained with a tolerance of 10 −3 Ryd/Bohr. For the bulk, unstrained structure, we took the experimental unit cell and optimized the positions of the atoms. An important structural parameter is the IrO 6 octahedral rotation angle θ, which we found to be 13.5 • as compared to the measured angle of 11.5 • .
For the epitaxially strained structure, the in-plane lattice constant was changed according to the strain parameter e xx ≡ (a − a bulk )/a bulk , while the out-of-plane lattice constant was determined from the condition that the cell volume is preserved under strain, so that e zz ≡ (c − c bulk )/c bulk = −2e xx . Structural parameters were then optimized for each strain condition, with the fixed lattice constants. We varied the in-plane strain e xx by ±2%, which is the same order of magnitude as the strains present in the experimental structures, e. g., e xx ≈ 2% for SIO grown on the GdScO 3 substrate and −2% for the NdGaO 3 substrate 8 . The computed structural parameters are shown in Table I and Fig. 2. The trend in the structural parameters shown here is similar to the results obtained from an earlier work using the dynamical mean-field theory 11 .
For test purposes, we also performed a second set of calculations, where we took the in-plane lattice constant to be the same as the substrate as usual, but used the measured c/a ratio from the experiments 8 . No significant changes in the bond angles and bond lengths were found.

III. DENSITY FUNCTIONAL BAND STRUCTURE
A. Electronic structure and Γ − X crossover The electron bands for the optimized structures under three different strain conditions are shown in Fig. 3, calculated within the local spin-density approximation with Coulomb interaction and spin-orbit coupling included (LSDA+SO+U). The Fermi surfaces for the doped SIO are shown in Fig. 4. The results show systematic changes of the band structure features in the gap region, some of which are shown in Fig. 5. The magnitude of the fundamental gap ∆ g increases with tensile strain, and so does the direct gap ∆ X at the X point, while the direct gap at the Γ point does not change very much. This trend in the gap values may be expected, since tensile strain would reduce the band widths of the LHB and the UHB, while the Coulomb U is relatively unchanged. The opposite happens for the compressive strain. The Γ point beats the trend because the valence top at Γ has a considerable contribution from the oxygen states. The approximately linear variation of the gap with strain is clearly visible in the partial density-of-states shown in Fig. 6. This trend of the gap variation with strain has already been noted in earlier works [11][12][13] .
Strain induced Γ -X crossover -An interesting feature of the band structure is the drastic change of the valence band maximum under strain, which would have a significant effect on the transport properties of the hole-doped structure. As can be seen from the band structure, Fig.  3, the conduction band minimum does not change under strain, always occurring at the M point, while it changes with strain. In the unstrained structure, DFT results show that the energy of the valence maximum at X and Γ are nearly the same. From our calculations, the Γ − X separation is E ΓX ≡ E Γ −E X = 0.01 eV, which compares well with the earlier DFT calculations in the literature, where E ΓX ≈ 0.01 − 0.07 eV 1,11,12,21 . Thus according to the DFT results, Γ is slightly below X, which is opposite to the ARPES measurement 22 , where the energy at X point is found to be above Γ by about 2 eV, i. e., E ΓX ≈ 0.2 eV 22 . For the tensile strain, the valence maximum at Γ increases in energy, while that of X goes down; For the compressive strain, the reverse happens, leading to a Γ − X crossover with strain.
The Γ − X crossover may be described within a tightbinding model involving the |J eff = 1/2 orbitals on the square lattice of Ir atoms and by including the dependence of the TB hopping integrals of the octahedral rotation angles which vary with strain. This TB model is developed in A. The TB model shows the correct trend with the result: ∆E ΓX (T B) ≈ 1.65 eV × e xx (Eq. A12) as compared to ∆E ΓX (DFT) ≈ 4.0 eV × e xx , obtained from the DFT results shown in Fig. 5. The TB analysis indicates that both the angle and distance changes with strain are important for the description of the Γ − X crossover. Based on the computed slope of E ΓX from DFT and the experimental Γ − X separation of 0.2 eV, we would need an estimated epitaxial compressive strain of e xx ≥ 3% in order to switch the valence band maximum from X to Γ. The Γ − X crossover is clearly seen from the Fermi surface plots shown in Fig. 4, where we have shown the Fermi surfaces for 5% dopant concentration.
For the case of electron doping, we find that the conduction band minimum always occurs at the M point in the Brillouin zone and the electron pocket is elliptical in shape.

B. Magnetic moment under strain
The calculated spin and orbital magnetic moments with and without strain are listed in Table II. We note that the total energy obtained with spin moments constrained along the planar direction (x) has a lower energy than the spin moments constrained alongẑ, consistent with thex orientation of the magnetic moments. Table  II shows the calculated magnetic moments using VASP with the LSDA+U+SO functional.
There has been much interest in the magnetic moments in SIO, in particular on the ratio µ l /µ s 11-13,23-25 , since a ratio of two indicates the spin-orbital entanglement of the wave function and a deviation from this value is indicative of a mixture between J eff = 1/2 and 3/2 sectors, as suggested from Table III. For the unstrained material (e xx = 0), the calculated total moment µ l + µ s = 0.38µ B is in good agreement with values obtained from the magnetic susceptibility measurements (0.5µ B ) 26 as well as from previous calculations (0.36µ B ) 1 . As already mentioned, the d 5 configuration of Ir can be thought of as a single hole in the t 2g manifold, leading to the cubic-field values, µ l = 2/3 µ B and µ s = 1/3 µ B , as seen from Table III. When a tetragonal field is present as in the case of SIO, the |J eff , m states get mixed among themselves, and the magnetic moments can be substantially altered from the cubic-field value.
From our calculations, we find µ l /µ s ≈ 2.6 (Table II) for the unstrained structure, in general agreement with earlier calculations 11,12 as well as with a recent measurement using non-resonant magnetic x-ray diffraction which obtained the value µ l /µ s ≈ 2.5. An earlier x-ray absorption measurement 23 yielded the ratio µ l /µ s ≈ 1.1; the reason for the discrepancy between the two measurements is unclear.
The magnetic moment ratio as well as its variation with strain can be approximately described by invoking a tetragonal crystal field for the single ion in the presence of the spin-orbit coupling. For more accurate description, a renormalized spin-orbit coupling has been invoked 11 . The results for the single ion in the tetragonal field are given in B. For the spin moment aligned alongx, which is the case for SIO, we have µ l /µ s = 2 + 4ξ/3, where ξ = ε/λ is the ratio of the tetragonal field to the spinorbit coupling strength λ ≈ 0.4 eV. For the unstrained . For compressive strain, the hole pocket is at the X point, while for tensile strain, the pocket shifts to Γ. For the unstrained case, valence top at Γ is only slightly below X, as indicated from the size of the two hole pockets in (b). The Fermi surface with 5% electron doping for the unstrained structure is shown in (d). It remains more or less unchanged with strain, unlike the hole case, with the elliptical electron pocket occurring at M . case, ε ≈ 0.14 eV 11,27 , so that the ratio µ l /µ s = 2.47, in reasonable agreement with the DFT result stated in Table II. If we use the results of Ref. 25 which suggest that the magnitude of ε increases (decreases) by about 0.1 eV for tensile (compressive) epitaxial strain of 2 %, then µ l /µ s = 2.8 for the tensile case and 2.1 for the compressive case, which more or less explains the DFT calculated trend shown in Table II.

IV. OPTICAL ABSORPTION
The change in the band structure with strain is reflected in the optical absorption spectrum. The basic quantity to compute is the imaginary part of the dimensionless dielectric constant ε 2 (ω), from which the real part ε 1 (ω) and the refractive index n(ω) are computed using the Kramers-Krönig relation. In the dipoleapproximation, the optical absorption coefficient α(ω) is given by where M cv (k) = ψ ck |p|ψ vk is the momentum matrix element between the conduction and the valence states (defined as unoccupied and occupied states, respectively), andê is the light polarization vector. A closely related subsidiary function, useful to the discussion of the optical absorption, is the joint density of states It is important to note that even though the JDOS(ω) scales up linearly with the size of the unit cell chosen, the optical absorption coefficient α(ω) and ε 2 (ω) are both independent of the size of the unit cell used in the band  Optical absorption coefficient α(ω) for bulk SIO without strain and for light polarization along the plane. The calculated spectra have been shifted to match with the measured optical gap. The origin of the three characteristic lowenergy peaks α, β, and γ is indicated in Fig. 8. calculation, as they must be. This is because a large number of transitions in the larger unit cell are simply forbidden. This is easily seen by realizing that if a larger unit cell is used in the calculation, then a large number of vertical transitions in the folded Brillouin zone (larger unit cell) is disallowed because the matrix element M cv (k) becomes zero from Bloch symmetry, even though the conduction and valence states both have nominally the same momentum k. Alternatively, one can start with the allowed optical transitions in the Brillouin zone corresponding to the smallest unit cell and then fold them into the smaller Brillouin zone corresponding to the larger unit cell and explicitly see that not all vertical transitions are allowed in the smaller Brillouin zone. The joint density of states JDOS (ω), in contrast, scales with the size of the unit cell, which can be seen by computing its integral over energy by taking advantage of the δ function in the definition. Note that the dipole approximation used in the expressions omits the local field and excitonic effects. Results are presented below for plane-polarized light with polarization direction in the plane or normal to the plane.

A. Unstrained bulk
The calculated absorption spectra for bulk SIO without strain are shown in Fig. (7) for light polarization along the plane (E plane), which is also compared to the measured data 8 . The spectra show three distinct peaks which can be understood from the transitions between the J eff = 1/2 and J eff = 3/2 states as shown in Fig. (8). The J eff = 3/2 states are split into two as seen from Table I due to the tetragonal crystal field V t = 2ε/3. In the experiment, only two distinct peaks are seen (these are labeled α and β), while the weaker γ peak is missing, possibly due to instrumental broadening.
We did not find any significant asymmetry in α(ω) for polarization directions within the plane. However, there is a large difference whether E plane or E z, as

Energy (eV)
Ir(5d)-t 2g spin ↓  The optical absorption in the low energy range occurs due to transitions within the Ir t 2g manifold, transitions that are however not dipole allowed due to the optical selection rules (δL = ±1 for plane-polarized light). In this subsection, we discuss the optical transition matrix elements and conclude that they become dipole allowed due to the admixture of the Ir p states into the Ir d bands. To get the dipole matrix elements, we need to estimate two things, as indicated from Eq. (4) below: (i) The amount of the Ir p -d admixture η, which we estimate from perturbation theory as well as from densityfunctional results and (ii) The optical matrix element for dipole-allowed transition between the Ir p and d orbitals p|ê·p|d , which we estimate from an effective hydrogenic model.
The d-d transitions are non-zero in the crystal only because there is deviation from the spherical symmetry, the so called crystal field, which also produces splitting between the t 2g and the e g states 28 . The crystal field not only splits the Ir (5d) orbitals, but it also mixes the Ir (5p) and (4f) orbitals into the 5d states, so that the d-d optical transitions have non-zero matrix element. Alternatively, such transitions may be equivalently described to be due to the transfer between Ir d atomic orbitals to the adjacent atoms because of the wave function overlap, e. g., within the linear combination of atomic orbitals (LCAO) model.
We now estimate the dipole matrix elements, which are non-zero due to the p − d and f − d mixing for the Ir atom. We first consider the p − d mixing, which as we shall see makes the larger contribution. Symbolically, this is given by The admixture η of the Ir p states into the Ir d bands, although small due to the large p − d energy difference, is nevertheless significant enough and is largely responsible for the optical transition in SIO. The magnitude of η can be estimated from the partial density of states (PDOS) in the band calculations in the energy range of the Ir d bands. With the definition for the PDOS, ρ d (ε) = i | ψ i |d | 2 δ(ε − ε i ) and ρ p (ε) = i | ψ i |p | 2 δ(ε−ε i ), and with the wave functions ψ ≈ d + ηp, one can estimate the average value of the admixture from the two PDOS in the energy range of interest. The PDOS, computed using the LAPW method, 29 are shown in Fig. 9 and we have obtained the value of η from the ratio of the integral of the two PDOS in the en-ergy window of the Ir t 2g bands, viz., from -2 to 1 eV. The result is: η 2 ≈ ρ p /ρ d ≈ 0.19 × 10 −2 or η ≈ 4%. We find a similar admixture of the 4f states, ρ f /ρ d ≈ 0.2 × 10 −2 , indicating a nearly equal 4f admixture η 4f ≈ 4% into the Ir (5d) bands.
Alternatively, the strength of the admixture η can be estimated from perturbation theory. From the second order perturbation theory, we have where the factor of three is from the degeneracy of the p states. If we further approximate p|V cf |d ∼ d|V cf |d ≈ ∆ cf , admittedly a crude approximation, and take the t 2g -e g splitting ∆ cf ∼ 2 eV, and use the result ε d − ε p ∼ 43 eV for Ir from the standard Atomic tables 30 , we get the perturbation theory result of η ≈ 14%, which is a factor of three too high as compared to the DFT result. However, given the crudeness of the approximation we made for the matrix elements in the perturbation theory, the order of magnitude agreement is reasonable. In the following, we shall use the DFT value for the p − d admixture η ≈ 4%. The second part is the estimation of the dipole matrix element M pd ≡ p|ê · p|d appearing in Eq. (4), which we do from an effective hydrogenic model 31,32 . In this model, the atomic wave functions are described as hydrogenic wave functions, but with an orbital dependent effective nuclear charge Z. Thus where ρ = (2Z/na 0 )r. Note that Z here is an effective atomic number, which takes into account the screening of the core electrons and depends on the principal quantum number n also. We take the value 31 Z = 18.7 appropriate for Ir n = 5 orbitals (5p and 5d). The integration can be performed analytically to yield the result M pd for the plane polarized light: in atomic units (viz., = 1, Bohr radius a 0 = 1, m e = 1/2, Energy unit = 1 Ryd). Note that for atoms, the direction of polarization doesn't matter due to spherical symmetry. Plugging in the estimated magnitudes of M pd and η in Eq. (4), we find the d − d matrix element due to the admixture of Ir (5p) into the Ir (5d) states A similar calculation for the matrix element due to the admixture with the Ir (4f) orbitals yields a number, which is an order of magnitude smaller than |M p dd |, viz., where we estimated the 4f admixture η f from Fig. (9) and used the effective atomic number Z = 38.3 for the Ir  4f orbitals 31 . Thus the Ir (4f) orbitals contribute a much smaller amount as compared to the Ir (5p) orbitals. The estimate given by Eq. (8), M dd ∼ 10 −2 , is very comparable to the average matrix element M cv ∼ 10 −2 obtained from the band calculations (Fig. 10), indicating that the admixture of the Ir (5p) orbitals is indeed responsible for the optical absorption.

C. Optical absorption under strain
Uniform strain -The calculated optical absorption for uniform compression (uniform pressure) is compared to the same for bulk SIO without any compression in Fig.  (11). Apart from a shift in the peak positions to lower energies, caused by a reduction in band gap due to the increase of the band width due to compression, we note that the overall absorption coefficient is increased for both light polarizations. This can be explained from our above argument of Ir 5p -5d admixture, because a larger crystal field upon compression mixes the Ir 5p more into the Ir 5d bands, making the optical matrix element larger. Indeed, following the same logic, as the lattice constant is increased, there is less and less admixture of the Ir (5p) orbitals, and in the limit of infinite lattice constant, the admixture vanishes (perfect spherical symmetry), so that optical absorption would be zero as would be expected for the d → d transition due to the dipole selection rules. Epitaxial strain -The change of the optical absorption with epitaxial strain is shown in Fig. (12). With compressive epitaxial strain (e xx < 0 and e zz > 0), the distances in the plane are reduced, leading to a larger optical absorption for E plane as argued for the uniform compression case. The same effect leads to an overall reduction of the optical absorption for E ẑ. For the absorption with E plane, all three peaks α, β, and γ are still there, but the peaks are red (blue) shifted with compressive (tensile) epitaxial strain, consistent with the band structure changes under strain shown in Fig. (3).
Optical anisotropy -The polarization dependence of the optical absorption has been reported in the literature 10 , where a significant anisotropy is found in the optical spectra. The absorption is significantly reduced for the polarization E ẑ as compared to E plane in the low energy region corresponding to the transition within the Ir t 2g manifold, while for the higher energy region, the difference in the relative strength is not as drastic. This is reproduced quite well from our calculations as seen from Fig. 13, where we also compare with the existing experiments. The calculated spectra is for the unstrained structure, while the experimental spectra is for a system with non-uniform strain, the only case for which experimental results are available for both polarizations. The anisotropy is directly attributable to the differences in the matrix element |ê · M cv (k)| (Eq. 2), since the variation of the refractive index n(ω) with energy or polarization is relatively weak, as indicated from Fig. 14.

V. SUMMARY
In summary, we studied the electronic properties and the optical absorption spectra of Sr 2 IrO 4 under epitaxial strain condition using density-functional methods. Systematic structural changes with strain were found including the staggered rotation angle θ, which was found to have important effect on the electronic structure. An interesting result is the Γ−X crossover of the J eff = 1/2 valence band maximum with strain, allowing for engineering of the hole pocket in the hole-doped material, with potentially drastic changes in the transport properties. A minimal tight-binding Hamiltonian was developed for the J eff = 1/2 sector, which is capable of qualitatively describing the important features of the band structure under strain, including the Γ − X crossover.
We calculated the optical spectra under strain, interpreting the results in terms of small admixtures of the Ir (5p) and (4f) states with the Ir (5d) bands, without which the optical transitions would be dipole forbidden. The calculated spectra were compared to the experiments, where available, and the observed anisotropy in the optical absorption, which is very strongly anisotropic for the low energy region but less strongly anisotropic for the high energy region (see Fig. 13), was correctly explained from the calculated results. Our work opens up the possibility of exploring the strain manipulation of the transport properties of epitaxially grown spin-orbit coupled Mott systems.
In this Appendix, we construct a minimal tight-binding model to describe the band structure within the J eff = 1/2 sector on the square lattice of Ir atoms, appropriate for SIO, in order to understand the Γ − X crossover of the valence band with epitaxial strain.
We are primarily interested in the J eff = 1/2 sector, viz., the lower and the upper Hubbard bands (LHB / UHB) that form the valence and the conduction bands in the gap region. A tight-binding description can follow two different paths: One is to keep the t 2g orbitals (in total six orbitals in the basis per Ir atom including spin) in the Hamiltonian and then fit the the two Hubbard bands in the DFT band structure and the other is to keep only the two J eff = 1/2 orbitals per Ir atom (the minimal model). For the fitting with the t 2g orbitals, the NN TB parameters that best fit the two Hubbard bands are (in eVs): V π = −0.20, V δ = 0.04, U = 0.78, λ = 0.4. In our discussions below, we employ the second approach, where we use the minimal TB model with just the two orbitals per Ir atom to describe the J eff = 1/2 sector.
We consider the Hubbard model on a square lattice with anti-ferromagnetic order as appropriate for SIO, keeping the two spin-orbital entangled J eff = 1/2 orbitals on each site, which we call e 1 and e 2 , defined with respect to the local octahedral axes, with the staggered rotations as indicated in Fig. (1), viz., The TB Hamiltonian is given by where c † iα creates an electron at the site i (which may be in sublattice A or B) and in the orbital e α , t ij is the hopping integral, which is non-zero only for hopping between the same type of orbitals e 1 or e 2 , U is the on-site Coulomb interaction, and the summation ij indicates sum over distinct pairs of bonds. Note that in general the hopping integral t ij are complex numbers as discussed below. We will retain hoppings up to the fourth NN as indicated in Fig. (15). Furthermore, we find that at least three (and sometimes four) NN hoppings need to be kept for an accurate description of the band structure in this minimal model for the J eff = 1/2 sector.
To obtain the hopping integrals between these orbitals, it is convenient to first obtain the integrals in the unrotated d basis from standard Tables 33 and then rotate the basis. Under rotation the angular momentum functions transform among one another, and in the present case, we have a site-dependent 5 × 5 rotation matrix appropriate for L = 2. We denote the unrotated basis as |α (xy, yz, zx, x 2 −y 2 , and 3z 2 −1, in that order), rotated basis by |α ′ on sublattice A and |α ′′ on sublattice B, and the corresponding rotation matrices that transform one basis into another by R ′ and R ′′ , i. e., |α ′ = R ′ |α and |α ′′ = R ′′ |α . The hopping integrals in the rotated basis are then given by H α ′ β ′′ ≡ α ′ |H|β ′′ = α|R ′T HR ′′ |β , or The rotation matrix R(θ) for L = 2 with rotation θ about theẑ axis is well known 34 Due to the staggered rotations, the hopping integrals between orbitals on the same sublattice and opposite sublattices are given using Eq. (A3) as follows: Denoting the two rotated J eff = 1/2 states corresponding to Eq. (A1) by | e 1A , | e 1B , | e 2A , and | e 2B , one can then obtain the hopping integrals between these set of orbitals from the standard d − d hopping integrals from Harrison's Tables 33 and using the above rotation matrices. With the standard notations for the direction cosines nlm for the distance vector joining the first atom to the second, these matrix elements are readily obtained. The hopping amplitude between the opposite sublattice is Notice that the hopping matrix is diagonal with complex elements for non-zero θ, which however can be made real by making a gauge transformation. For hopping between the same sublattice, the hopping matrix is diagonal, but with real elements this time: e iA |H| e jA = tδ ij , FIG. 15: Hopping integrals between the | e1 orbitals, which can be obtained from Eqs. (A5) and (A6). The samesublattice hoppings t2 and t3 are θ independent and real, while the opposite-sublattice hoppings, t1 and t4, depend on θ and are complex if θ = 0. Hopping integrals for | e2 orbitals are complex conjugate of those for the | e1 orbitals.
We make a further simplification by taking V σ = −(3/2)V π following Harrison's scaling 33 in order to construct a minimal model. The same expression Eq. (A6) is valid for hopping between B sublattice atoms e iB |H| e jB , except that the sign of θ is changed due to the staggered rotations. Sometimes, in the literature, a simplistic angle dependence for the hopping h(θ) = h 0 cos θ is used. However, as the two equations above demonstrate, the angle dependence is more complicated and cannot be written down as cos θ even to the lowest order in the angle. In our TB model, we have retained up to four NN hoppings, as indicated in Fig. 15.
Gauge transformation -As seen from Eq. (A5), the hopping integrals are in general complex, and of course there is no problem working with the complex hoppings. However, it is convenient to make the NN integral real by a gauge transformation 35 . In the transformation, one simply adds a multiplicative phase factor to the orbital definitions, viz., where ε j = ±1 for j = A/B sublattices. With the choice φ = tan −1 (t ′ /t), where t, t ′ are given in Eq. (A5) for the 1NN hopping, the new hopping integral becomes real, viz., t + it ′ →t, wheret = (t 2 + t ′2 ) 1/2 . This transformation leaves all same-sublattice hoppings real, while the opposite-sublattice hoppings including t 4 and beyond continue to remain complex. One can choose to work with such complex hoppings, or else, simply ignore the imaginary parts, since hoppings between far-away neighbors are small anyway. We denote the gauge-transformed basis (A7) as |e 1 and |e 2 in Eq. (A1) and the hopping amplitudes (real) as t i as indicated in Fig. (15). The antiferromagnetic lattice structure together with the fact that | e 1 and | e 2 subspaces do not mix leads to the 2 × 2 TB Hamiltonian in the momentum space where h 11 = 4t 2 cos k x cos k y + 2t 3 (cos 2k x + cos 2k y ), h 12 = 2t 1 (cos k x + cos k y ) + 4t 4 (cos 2k x cos k y + cos 2k y cos k x ), a 0 = 1 here, and ∆ = U/2 is the staggered field. Diagonalization readily yields the energies of the upper and lower J eff = 1/2 Hubbard bands to be The TB parameters obtained from fitting Eq. (A9) to the density-functional bands for the optimized structure with no strain are, in units of eV: U = 0.65, t 1 = −0.095, t 2 = 0.015, t 3 = 0.035, and t 4 = 0.01. The TB fit to the DFT bands are shown in Fig. (3) as dotted lines. An important feature of the band structure is the occurrence of the conduction minimum at the M point of the Brillouin zone, which the electrons would occupy in the electron doped system. It can be easily shown that for this to happen, the TB parameters must satisfy the condition 36 which is clearly satisfied by our parameters given above. Furthermore, the eccentricity of the elliptical energy contours around the M point in the conduction band, as seen in Fig. 4, is given by 37 We determine the TB parameters under strain condition with the following ansatz. Using the calculated angle θ ≈ 13 • and the hopping expressions Eqs. (A5) and (A6), there is a one to one correspondence between t i and V i π for the i-th neighbor. We obtain the values of V i π for different neighbors under no strain condition. In order to compute the hopping integrals under strain conditions, we back substitute V i π and the DFT optimized angles under strain into Eqs. (A5) and (A6). In addition to the change of the angles, distances between atoms also change. Taking the variation of V π to follow Harrison's R −5 scaling with distance R and including the effect of θ from Eqs. (A5) and (A6), we can compute the TB hopping parameters t i under the epitaxial strain conditions. We assume that strain enters the hopping parameters only via change of the hopping distances and the rotation angle θ.
With this simple TB model, we find that the important trends of the band structure under strain are described qualitatively correctly, though not quantitatively with a factor of two to three discrepancy as compared to the DFT results.
In particular, the model can qualitatively describe the Γ − X crossover of the valence band top under strain (hole pocket in the doped structure), predicted from the DFT calculations. From the TB energy expression Eq. (A9), we readily find the Γ − X energy difference for the valence band (A11) For our TB parameters, this is almost zero consistent with the DFT results. The strain dependence of E ΓX is obtained from the TB expression (A11), viz., where we have used the fact that the distance change between the atoms in the plane is simply ∆R/R = e xx , and also the DFT result for the change in angle with strain, viz., ∆θ ≈ (−1.5 rad) e xx , obtained from Table  I assuming a linear strain dependence, since the strain is small. It turns out that the two terms in Eq. (A12) contribute nearly half each to the final result, so that both the angle and distance changes are important for the description of the Γ − X crossover.
Appendix B: Single ion in a tetragonal field with spin-orbit coupling It is instructive to examine the eigenstates of the single site Hamiltonian for the d orbitals in the presence of the SOC and a tetragonal crystal field, which this Appendix deals with.
We assume a large cubic crystal field splitting ∆ cf → ∞, so that the e g and the t 2g sectors don't mix. Within the t 2g sector, the Hamiltonian is given by By diagonalizing the Hamiltonian and keeping the tetragonal field to the lowest order, we obtain the eigen-  III: Energies E and wave functions of the Ir d atomic states in the presence of spin-orbit coupling λ L · S and cubic (∆ cf ) as well as tetragonal crystal fields (ε, δ assumed to be ≪ λ and treated perturbatively, keeping only the linear terms in ǫ/λ). The standard |J eff , m labels for the t2g states are also indicated along with the expectation values Lz and 2Sz , if the net spin is alongẑ, and Lx and 2Sx , if the net spin is alongx. The quantity ξ ≡ ǫ/λ is the ratio of the tetragonal field ǫ to the spin-orbit coupling constant λ.
Cubic Field (O h ) Tetragonal field (D 4h ) orbital E |J eff , m E wave functions net spin alongẑ net spin alongx 2Sz Lz values and eigenfuctions, which are summarized in Table  III including some relevant expectation values. Note that the wave functions |J eff , m listed in Table III are not necessarily eigenstates of J 2 and J z , because they were not obtained by diagonalizing the full 10 × 10 λ L · S matrix, but rather just the 6 × 6 matrix in the t 2g sector. The basic electronic structure of the iridates is determined by the Ir 4+ ions with the 5d 5 configuration and a large crystal field due to the oxygen octahedra. The crystal field combined with the SOC results in the spinorbital entangled states as summarized in Table III. The cubic crystal field splits the 5d states into e g plus t 2g states. With the SOC included, the six-fold degenerate t 2g states (including spin) split into a two-fold J eff = 1/2 (higher energy) and a four-fold J eff = 3/2 (lower energy) state. With the 5d 5 configuration, four electrons fill the J eff = 3/2 state, while the remaining electron occupies the lower Hubbard band separated from the upper Hubbard band by the Coulomb interaction U thus resulting in a SOC-induced Mott insulator. The tetragonal field has been modeled in Table III by asymmetric on-site energies for the Ir d orbitals, viz., ε = ε(xz/yz) − ε(xy) and δ = ε(3z 2 − 1) − ε(x 2 − y 2 ). The d 5 configuration is therefore equivalent to a single hole t 1 2g with the |1/2, 1/2 configuration.
In Table III, we have also listed the angular momentum expectation values, which will be useful for the interpretation of the computed magnetic moments discussed in Section III B. The Table lists L z and 2S z for the spin state alongẑ; It is straightforward to show that L x = L y = S x = S y = 0 in this case. Also listed are the expectation values L x and 2S x for the spin state alonĝ x direction (which happens to be the case for the ground state of SIO); Again, L y = L z = S y = S z = 0 for this case. One point to note here is that the tetragonal field ξ mixes up the J eff = 1/2 and 3/2 wave functions. For spin magnetization alongẑ, the expectation values of the magnetic moments are µ z l = µ B ψ 1/2 |L z |ψ 1/2 = 2/3 + 8ξ/27, (B2) µ z s = µ B ψ 1/2 |2S z |ψ 1/2 = 1/3 + 16ξ/27, where |ψ 1/2 = |1/2, 1/2 , the Bohr magneton µ B = 1, and other components of the magnetic moments are zero, leading to the ratio of the magnetic moments to be µ z l /µ z s = 2 − 8ξ/3.