Strain-induced superconductivity in Sr₂IrO₄

Lena Engström [email protected] Department of Physics and the Centre for the Physics of Materials, McGill University, Montréal, Québec, H3A 2T8, Canada Département de Physique, Université de Montréal, Montréal, Québec, H3C 3J7, Canada Regroupement Québécois sur les Matériaux de Pointe (RQMP) Chia-Chuan Liu Département de Physique, Université de Montréal, Montréal, Québec, H3C 3J7, Canada William Witczak-Krempa Département de Physique, Université de Montréal, Montréal, Québec, H3C 3J7, Canada Regroupement Québécois sur les Matériaux de Pointe (RQMP) Centre de Recherches Mathématiques, Université de Montréal, Montréal, Québec, HC3 3J7, Canada Institut Courtois, Université de Montréal, Montréal, Québec, H2V 0B3, Canada T. Pereg-Barnea Department of Physics and the Centre for the Physics of Materials, McGill University, Montréal, Québec, H3A 2T8, Canada Regroupement Québécois sur les Matériaux de Pointe (RQMP)

Abstract

Multi-orbital quantum materials with strong interactions can host a variety of novel phases. In this work we study the possibility of interaction-driven superconductivity in the iridate compound Sr₂IrO₄ under strain and doping. We find numerous regimes of strain-induced superconductivity in which the pairing structure depends on model parameters. Spin-fluctuation mediated superconductivity is modeled by a Hubbard-Kanamori model with an effective particle-particle interaction, calculated via the random phase approximation. Magnetic orders are found using the Stoner criterion. The most likely superconducting order we find has $d$ -wave pairing, predominantly in the total angular momentum, $J=1/2$ states. Moreover, an $s_{\pm}$ -order which mixes different bands is found at high Hund’s coupling, and at high strain anisotropic $s$ - and $d$ -wave orders emerge. Finally, we show that in a fine-tuned region of parameters a spin-triplet $p$ -wave order exists. The combination of strong spin-orbit coupling, interactions, and a sensitivity of the band structure to strain proves a fruitful avenue for engineering new quantum phases.

I Introduction

The iridates display a rich phase diagram due to an interplay between strong correlations, spin-orbit coupling, and crystal field effects, all acting on multiple $d$ -orbitals which in some cases lead to multiple Fermi surfaces. Various iridates show novel phenomena such as Kitaev and Weyl physics [1, 2, 3]. The first compound in the family of Ruddlesden-Popper perovskite strontium iridates, Sr₂IrO₄, consists of stacked quasi-2d layers. In each layer, the iridium atoms form a square lattice, and are surrounded by octahedra of oxygen atoms. As shown in Fig. 1, every other IrO₆ octahedron is rotated by an angle of $\phi=\pm\phi_{\epsilon=0}\approx\pm 12^{\circ}$ relative to the iridate lattice and we therefore use a two-site basis [4]. In Sr₂IrO₄ the three $t_{2g}$ orbitals are located close to the Fermi surface. However, due to a strong spin-orbit coupling the resulting band structure is better characterized by the on-site total angular momentum $J$ eigenstates, with $J=1/2$ states having $J_{z}=\pm 1/2$ and $J=3/2$ having the projections along the $z$ -axis $J_{z}=\pm 1/2,\pm 3/2$ . We refer to this basis as the $j$ -states. For the undoped compound the electron filling is $n=5$ , out of the 6 $t_{2g}$ bands per site, and an antiferromagnetic insulator is found up to $T_{\textrm{N}}=240$ K [5]. For this state, the two bands of mainly $j=1/2$ character are half-filled and located close to the Fermi surface while additional bands of $j=3/2$ character are further from the Fermi level [6].

Refer to caption — Figure 1: Single plane of Sr₂IrO₄ where the iridium atoms (gray) form a square lattice, and are surrounded by oxygen octahedra. There are two sites in the unit cell $s=A,B$ with staggered rotations $\phi_{\epsilon}$ of the octahedra, resulting in a bond angle $\theta_{\epsilon}=180^{\circ}-2\phi_{\epsilon}$ . The staggered rotation angle increases with compressive strain $\phi_{\epsilon}>\phi_{\epsilon=0}$ , with $\phi_{\epsilon=0}\approx 12^{\circ}$ .

Experimentally, Sr₂IrO₄ shows Fermi arcs and a pseudogap under electron doping [7, 8] as well as non-Fermi liquid behavior under hole doping [8]. Superconductivity has been predicted for both types of charge doping of Sr₂IrO₄ in multiple theoretical works [9, 10, 11, 12, 13, 14]. At a first approximation the band structure and the interactions suggest a direct parallel to known high $T_{c}$ superconductors, the cuprates. Indeed, projecting the Hamiltonian on the $j=1/2$ states results in a model similar to those used to describe the cuprates [15] and a simplified one-band model therefore predicts $d$ -wave superconductivity for electron doped Sr₂IrO₄ [9]. However, theoretical studies of Sr₂IrO₄ found that $d$ -wave superconductivity is likely to arise only in a limited range of interaction parameters. At hole doping, additional pockets of $j=3/2$ character appear at the Fermi level. In this region of the phase diagram, previous studies have predicted Sr₂IrO₄ to have either multi-band $s_{\pm}$ -wave or a $p$ -wave pairing. However, as of yet no superconducting order has been experimentally confirmed for Sr₂IrO₄ when chemical doping is the only tuneable parameter [16, 17, 18, 19]. The question then remains whether there could be a tuning parameter that would make superconductivity more favorable.

The local staggered rotations of iridium sites introduce an in-plane translation symmetry breaking accompanied by additional hybridization of orbitals. These effects have been previously ignored in multi-orbital models of iridate superconductivity [20]. The rotations increase under compression and as a result the hopping between orbitals at neighboring sites is modified to reflect the new geometry [21, 4, 22, 23, 24]. Moreover, the orbitals are modified by different amounts such that the bands belonging to the $j=3/2$ subspace move closer to the Fermi surface [25]. Naturally, the number of Fermi pockets and their orbital composition are important factors for superconductivity. As the undoped Sr₂IrO₄ is an antiferromagnetic insulator, a prerequisite for any superconducting order is that it must exist in a regime where the system is no longer magnetic. Several experimental studies have shown that by growing Sr₂IrO₄ on a substrate with mismatched lattice parameters, the induced compressive epitaxial strain significantly suppresses the magnetic order [26, 27, 28, 29, 30, 31]. In a variety of known superconductors biaxial, either compressive or tensile, strain has proven to increase the critical temperature [32, 33, 34, 35, 36] or to induce a superconductivity/insulator phase transition [37, 38]. In the current work compressive strain is suggested to induce the same type of phase transition and could thus expand the region of doping where superconductivity can be observed. The purpose of our study is therefore to determine if superconductivity is more likely when strain is applied. More precisely, we aim to answer the following two questions. First, are there regimes of applied strain where a superconducting order is possible? Second, does the in-plane symmetry breaking due to the rotations result in different superconducting orders?

Experiments in undoped Sr₂IrO₄ under high hydrostatic pressure have been performed in recent years [39, 40, 41]. While both hydrostatic pressure and epitaxial strain change the interatomic distances, they do so in different ways. The experiments approximate the pressure to strain conversion to be $(\Delta a/a)/\Delta P=-0.146\%/$ GPa [39]. A transition into a non-magnetic insulating state occurs under a compression of 17GPa. At sufficient pressure beyond that the resistivity shows a rapid increase accompanied by a pressure-induced structural phase transition [39, 40]. However, it is important to remember that while the effect on in-plane distances is similar, hydrostatic pressure decreases the inter-layer distance while epitaxial strain increases it [31]. In the hydrostatic pressure, the $c$ -axis compression increases interactions between perovskite layers while the epitaxial strain does not. A persistent insulating state is thus not expected in the realistic regimes of our phase diagrams. In addition, our region of interest is for charge doping, where the insulating nature of the compound is weaker.

Spin-fluctuations are believed to be able to mediate superconductivity in the iridates [12, 14]. Multi-orbital superconductivity has successfully been modeled with spin-fluctuations in other families of materials such as ruthenates [Romer2019, Kaser2022, Gingras2022] and iron-based superconductors [Ikeda2010, Schrodi2020]. In this work, a linearized superconducting gap equation (in the static limit) is solved to find regimes where superconductivity is possible. We consider a multi-orbital Hubbard-Kanamori model of Sr₂IrO₄ in a rotated two-site basis. The spin susceptibility is calculated via the random phase approximation (RPA). The spin fluctuations are thus dependent on the staggered sublattice rotations. As the rotations increase with an increasing strain and the RPA susceptibility is used to derive the effective particle-particle interaction, the interaction is dependent on the strain. This in turn results in a strain-dependent linearized gap equation for the superconducting order. Magnetic orders can be identified via the RPA susceptibility. We find a large region of strain-induced superconductivity, as well as several possible magnetic orders. The different types of superconducting order are either mediated by spin or pseudospin fluctuations. We find that as the compressive strain is increased the fluctuations become more spin-like in character. Although several types of fluctuations compete in parts of the calculated phase diagrams, the most prevalent type is antiferromagnetic fluctuations in the $j=1/2$ state. These pseudospin fluctuations can mediate a $d$ -wave order. Longer range fluctuations in the spin basis instead mediate the $s_{\pm}$ -wave superconductivity. For high compressive strain, intraorbital spin fluctuations can become large enough to mediate anisotropic superconducting orders. All superconductivity in the calculated phase diagrams are mediated by fluctuations of spins oriented in-plane. However, there exists regions where ferromagnetic out-of-plane fluctuations are of equal size. The ferromagnetic fluctuations are found to mediate an odd parity $p$ -wave order.

The paper is structured as follows. In section II, we introduce the underlying tight-binding Hamiltonian of Sr₂IrO₄ and how the compressive strain is modeled. In subsection II.2, the model for superconductivity mediated by spin-fluctuations is introduced. The resulting phase diagram is presented in section III, for realistic values of model parameters. Additional phase diagrams are shown for a wider variety of possible values for the Hund’s and spin-orbit coupling. We then analyze the nature of the magnetic fluctuations in section IV. The types of superconducting orders, and the fluctuations believed to mediate them, are detailed in section V. Finally, in section VI we discuss the experimental possibilities of strain-induced superconductivity, and signatures of the found orders.

II Model and methods

II.1 Kinetic Hamiltonian with rotations

The band structure of Sr₂IrO₄ can be modeled by a spin-orbit coupled tight-binding Hamiltonian:

H=H_{\textrm{kin}}+H_{\textrm{SOC}}

(1)

We consider a 2-site orbital-spin basis: $\boldsymbol{c}=$ $(c_{\boldsymbol{k},A,yz,\uparrow},$ $c_{\boldsymbol{k},A,yz,\downarrow},$ $c_{\boldsymbol{k},A,xz,\uparrow},$ $c_{\boldsymbol{k},A,xz,\downarrow},$ $c_{\boldsymbol{k},A,xy,\uparrow},$ $c_{\boldsymbol{k},A,xy,\downarrow},$ $c_{\boldsymbol{k},B,yz,\uparrow},$ $c_{\boldsymbol{k},B,yz,\downarrow},$ $c_{\boldsymbol{k},B,xz,\uparrow},$ $c_{\boldsymbol{k},B,xz,\downarrow},$ $c_{\boldsymbol{k},B,xy,\uparrow},$ $c_{\boldsymbol{k},B,xy,\downarrow})$ . For each spin $\sigma=\uparrow,\downarrow$ the kinetic terms have intra- and inter-sublattice hopping:

H_{\textrm{kin}}=\left(\begin{array}[]{c c }H_{AA}&e^{ik_{x}}H_{AB}\\ e^{-ik_{x}}H_{AB}^{\dagger}&H_{BB}\end{array}\right)

(2)

H_{AA}=\left(\begin{array}[]{c c c}\epsilon_{d}&\epsilon_{1d}&0\\ \epsilon_{1d}&\epsilon_{d}&0\\ 0&0&\epsilon^{xy}_{d}\end{array}\right),H_{AB}=\left(\begin{array}[]{c c c}\epsilon_{yz}&-\epsilon_{rot}&0\\ \epsilon_{rot}&\epsilon_{xz}&0\\ 0&0&\epsilon_{xy}\end{array}\right)

(3)

and $H_{BB}=H_{AA}$ . The factor $e^{ik_{x}}$ arises from the choice of unit cell, where the two sublattice sites are chosen as in Fig. 1 and the lattice spacing, $a$ , is set to $1$ . The hopping terms are

$\displaystyle\epsilon_{xy}$	$\displaystyle=2t\left(\cos k_{x}+\cos k_{y}\right)$	(4)
$\displaystyle\epsilon_{yz}$	$\displaystyle=2\left(t_{\delta}\cos k_{x}+t_{1}\cos k_{y}\right)$
$\displaystyle\epsilon_{xz}$	$\displaystyle=2\left(t_{1}\cos k_{x}+t_{\delta}\cos k_{y}\right)$
$\displaystyle\epsilon_{rot}$	$\displaystyle=2t^{\prime}\left(\cos k_{x}+\cos k_{y}\right)$
$\displaystyle\epsilon^{xy}_{d}$	$\displaystyle=4t_{n}\cos k_{x}\cos k_{y}+\mu_{xy}$
$\displaystyle\epsilon_{1d}$	$\displaystyle=4t_{1d}\sin k_{x}\sin k_{y}$
$\displaystyle\epsilon_{d}$	$\displaystyle=4t_{nd}\cos k_{x}\cos k_{y}.$

The hopping values are $(t,t_{1},t_{\delta},t^{\prime},$ $t_{n},t_{1d},t_{nd},\mu_{xy})=$ $-0.36(1,0.882,0.260,0.199,0.559,$ $0.019,-0.010,0.7)$ eV, when no strain is applied. The staggered rotations result in the non-zero inter-orbital hopping, $\epsilon_{rot}$ and $\epsilon_{1d}$ , between the $yz$ - and $xz$ -orbitals. When the compressive strain, $\epsilon<0$ , is increased the hopping parameters change. We use a linear strain dependence as in Ref. [25], following data by Ref. [31]:

$\displaystyle t(\epsilon)$	$\displaystyle=t\left(1+\rho\epsilon\right)$	(5)
$\displaystyle t_{1}(\epsilon)$	$\displaystyle=t_{1}\left(1+\rho_{1}\epsilon\right)$
$\displaystyle t^{\prime}(\epsilon)$	$\displaystyle=t^{\prime}\left(1+\rho^{\prime}\epsilon\right)$
$\displaystyle t_{n}(\epsilon)$	$\displaystyle=t_{n}\left(1+\rho_{n}\epsilon\right)$
$\displaystyle t_{\delta}(\epsilon)$	$\displaystyle=t_{\delta}\left(1+\rho_{\delta}\epsilon\right)$
$\displaystyle t_{1d}(\epsilon)$	$\displaystyle=t_{1d}\left(1+\rho_{1d}\epsilon\right)$
$\displaystyle t_{nd}(\epsilon)$	$\displaystyle=t_{nd}\left(1+\rho_{nd}\epsilon\right)$
$\displaystyle\phi_{\epsilon}$	$\displaystyle=\phi_{\epsilon=0}\left(1+\rho_{\phi}\epsilon\right)$

with $(\rho,\rho_{1},\rho^{\prime},\rho_{n},\rho_{\delta},\rho_{1d},\rho_{nd},\rho_{\phi})=$ (0.014,-0.251, -0.309, -0.048, 0, 0,-0.02,-0.085). $\phi_{\epsilon}$ is angle of the staggered rotations. The approximation considers not only rigid rotations of the oxygen octahedra but also changes in bond lengths, as was recently proposed to be a more accurate description of the strain effect in Ref. [Parschke2022]. The strain, and the associated rotations of the octahedra, increase the inter-orbital $t^{\prime}$ and intra-orbital $t_{1}$ . On the other hand hopping within the $xy$ -orbital decreases. The effect of strain on the tetragonal splitting $\mu_{xy}$ is not discussed here and is deferred to Appendix C.

It has been well-established that the spin-orbit coupling in Sr₂IrO₄ is large enough for each band to have a clear character of either total angular momenta $j=1/2$ or $j=3/2$ . The atomic SOC is

H_{\textrm{SOC}}=\frac{\lambda}{2}\sum_{\alpha\beta,\sigma\sigma^{\prime}}\sum_{s=A,B}\boldsymbol{L}_{\alpha\beta}\cdot\boldsymbol{\sigma}_{\sigma\sigma^{\prime}}c^{\dagger}_{\boldsymbol{k}s\alpha\sigma}c_{\boldsymbol{k}s\beta\sigma^{\prime}}

(6)

where $\boldsymbol{\sigma}=\left(\sigma^{x},\sigma^{y},\sigma^{z}\right)$ are the Pauli matrices in the spin basis $\sigma=\uparrow,\downarrow$ with respect to the $z$ -direction, and

\boldsymbol{L}=\left(\begin{bmatrix}0&0&0\\ 0&0&-i\\ 0&i&0\end{bmatrix},\begin{bmatrix}0&0&i\\ 0&0&0\\ -i&0&0\end{bmatrix},\begin{bmatrix}0&-i&0\\ i&0&0\\ 0&0&0\end{bmatrix}\right).

(7)

The eigenstates of $H_{\textrm{SOC}}$ are the $j$ -states, with associated annihilation operators: $\boldsymbol{a}_{m,\tau}=\sum_{\alpha,\sigma}M_{(m,\tau),(\alpha,\sigma)}\boldsymbol{c}_{\alpha,\sigma}$ . Where

M=\left(\begin{array}[]{cccccc}0&\frac{1}{\sqrt{3}}&0&-\frac{i}{\sqrt{3}}&\frac{1}{\sqrt{3}}&0\\ \frac{1}{\sqrt{3}}&0&\frac{i}{\sqrt{3}}&0&0&-\frac{1}{\sqrt{3}}\\ 0&\frac{1}{\sqrt{6}}&0&-\frac{i}{\sqrt{6}}&-\sqrt{\frac{2}{3}}&0\\ -\frac{1}{\sqrt{6}}&0&-\frac{i}{\sqrt{6}}&0&0&-\sqrt{\frac{2}{3}}\\ 0&-\frac{1}{\sqrt{2}}&0&-\frac{i}{\sqrt{2}}&0&0\\ \frac{1}{\sqrt{2}}&0&-\frac{i}{\sqrt{2}}&0&0&0\end{array}\right)

(8)

with $\boldsymbol{c}_{\alpha,\sigma}=(c_{yz,\uparrow},c_{yz,\downarrow},c_{xz,\uparrow},c_{xz,\downarrow},c_{xy,\uparrow},c_{xy,\downarrow})$ and $\boldsymbol{a}_{m,\tau}=(a_{1,+},a_{1,-},a_{2,+},a_{2,-},a_{3,+},a_{3,-})$ . In this basis each site has the states $a_{m,\tau}$ , where $m$ denotes the total angular momentum and its $z$ -axis projection $\left(j,j_{z}\right)$ such that $1=\left(\frac{1}{2},\pm\frac{1}{2}\right),2=\left(\frac{3}{2},\pm\frac{1}{2}\right),3=\left(\frac{3}{2},\pm\frac{3}{2}\right)$ and the projections along the $z$ -axis are labeled by $\tau=\pm$ . The projections $\tau$ can be treated as pseudospins, here not mixing sublattice and spin degrees of freedom but orbital and spin [Sigrist1991]. The total of 12 bands $b_{\boldsymbol{k},n}$ have eigenvalues $\xi_{\boldsymbol{k},n}$ and are connected to the orbital basis via $c_{\boldsymbol{k},j}=\sum_{n}U_{\boldsymbol{k},jn}b_{\boldsymbol{k},n}$ . The only spin-mixing in the Hamiltonian comes from the atomic spin-orbit coupling and all hopping terms are pseudospin-conserving. The non-interacting Hamiltonian is therefore separable into pseudospin $\tau=+,-$ sectors, containing the states $\{(yz,\downarrow),(xz,\downarrow),(xy,\uparrow)\}$ and $\{(yz,\uparrow),(xz,\uparrow),(xy,\downarrow)\}$ respectively. Therefore the 12 bands can be described by 6 bands in each pseudospin-sector $n_{\tau}$ .

As the strain modifies the hopping parameters, both the shape and the number of pockets at the Fermi surface changes. Therefore the Fermi surface is different at every point in the phase diagrams. As can be seen in Figs. 2, there are pockets belonging to the bands of $\left(j,j_{z}\right)=\left(\frac{1}{2},\pm\frac{1}{2}\right)$ present at every point in the phase diagram. For fillings corresponding to hole doping, $n<5$ , another pocket with $\left(j,j_{z}\right)=\left(\frac{3}{2},\pm\frac{3}{2}\right)$ appears around $(k_{1},k_{2})=(0,0)$ . The strain increases the bandwidth and cause the $(\frac{3}{2},\pm\frac{3}{2})$ pocket to appear for all doping values. As shown further in Appendix B the size of the $j=1/2$ electron pocket increases with strain.

II.2 Spin-fluctuation mediated superconductivity

To analyze the possibility of spin fluctuation mediated superconductivity, we solve a linearized gap equation in the static limit and normal state. A similar calculation has been performed previously for Sr₂IrO₄ in Ref. [14], for a model without staggered rotations or strain. In general, the particle-hole and particle-particle self-energies are defined by the Dyson-Gorkov equation:

\boldsymbol{G}(k)=\boldsymbol{G}^{(0)}(k)+\boldsymbol{G}^{(0)}(k)\boldsymbol{\Sigma}(k)\boldsymbol{G}(k)

(9)

with

\displaystyle\boldsymbol{G}(k)=\left[\begin{array}[]{c c}G(k)&F(k)\\ \bar{F}(k)&\bar{G}(k)\end{array}\right],\qquad\boldsymbol{\Sigma}(k)=\left[\begin{array}[]{c c}\Sigma(k)&\Delta(k)\\ \bar{\Delta}(k)&\bar{\Sigma}(k)\end{array}\right]

(10)

\displaystyle\boldsymbol{G}^{(0)}(k)=\left[\begin{array}[]{c c}G^{(0)}(k)&0\\ 0&\bar{G}^{(0)}(k)\end{array}\right]

(11)

where $G(k)$ is the particle-hole and $\bar{G}(k)$ the hole-particle Green’s functions

G_{ab}(k)=-\int_{0}^{\beta}d\tau e^{i\omega_{n}\tau}\langle T_{\tau}c_{\boldsymbol{k}a}(\tau)c_{\boldsymbol{k}b}^{\dagger}(0)\rangle

(12)

\bar{G}_{ab}(k)=-\int_{0}^{\beta}d\tau e^{i\omega_{n}\tau}\langle T_{\tau}c_{-\boldsymbol{k}a}^{\dagger}(\tau)c_{-\boldsymbol{k}b}(0)\rangle

(13)

for crystal momentum $\boldsymbol{k}$ and fermionic Matsubara frequencies $\omega_{n}=(2n+1)\pi/\beta$ , with $\beta=\hbar/(k_{\textrm{B}}T)$ . In the orbital-spin-sublattice basis $a=(\alpha,\sigma,s)$ with orbital $\alpha$ , spin $\sigma$ , and sublattice $s$ . The non-interacting Green’s functions are calculated with the operators $c_{\boldsymbol{k}a}$ for the Hamiltonian in Eq. (1). In the band ( $n$ ) and orbital-spin-sublattice basis respectively:

G^{(0)}_{n}(\boldsymbol{k},i\omega_{n})=\left[i\omega_{n}-\xi_{\boldsymbol{k},n}\right]^{-1}

(14)

G^{(0)}_{ab}(k)=\sum_{n}U_{\boldsymbol{k},an}G^{(0)}_{n}(k)U^{\dagger}_{\boldsymbol{k},nb}.

(15)

For the non-interacting model $G^{(0)}_{ab}(k)=-\bar{G}^{(0)}_{ba}(-k)$ . Under the approximations of this work, only the particle-particle self-energy is calculated. In a multi-orbital system, the non-interacting particle-hole susceptibility $\hat{\chi}_{0}(k)=\hat{\chi}_{0}^{ph}(k)$ is a tensor calculated from the non-interacting Green’s functions in the $N_{k}\times N_{k}$ -lattice we are considering:

\chi^{ph}_{0,abcd}(q)=\frac{1}{N_{k}^{2}\beta}\sum_{k}G^{(0)}_{ac}(q+k)\bar{G}^{(0)}_{bd}(-k)

(16)

Using the well-known summation for the Lindhard function over all fermionic Matsubara frequencies $\omega_{n}$ and the analytic continuation:

		$\displaystyle\chi_{0,abcd}(\boldsymbol{q},i\omega_{n}\rightarrow i\delta)=$		(17)
		$\displaystyle\frac{1}{N_{k}^{2}}\sum_{\boldsymbol{k},n,n^{\prime}}\left[U_{\boldsymbol{k}+\boldsymbol{q}}\right]_{an}\left[U^{\dagger}_{\boldsymbol{k}+\boldsymbol{q}}\right]_{nc}\left[U_{\boldsymbol{k}}\right]_{dn^{\prime}}\left[U^{\dagger}_{\boldsymbol{k}}\right]_{n^{\prime}b}$
		$\displaystyle\times\frac{f(\xi_{\boldsymbol{k}+\boldsymbol{q},n},T)-f(\xi_{\boldsymbol{k},n^{\prime}},T)}{i\delta-\left(\xi_{\boldsymbol{k}+\boldsymbol{q},n}-\xi_{\boldsymbol{k},n^{\prime}}\right)}$

where $\delta$ is a small number, here set to $10^{-4}$ eV. $f(\xi_{\boldsymbol{k},n},T)$ is the Fermi-Dirac distribution at temperature $T$ . The tensor $\hat{\chi}_{0}(q)$ can be written as a rank 4 tensor. The number of spins $N_{\sigma}=2$ , orbitals $N_{o}=3$ , and sublattices $N_{s}=2$ , results in a tensor of dimension $12\times 12\times 12\times 12$ . However, to separate the spin degrees of freedom we can reshape the tensor to have the dimension $N_{\sigma}^{2}\times N_{\sigma}^{2}\times(N_{o}N_{s})^{2}\times(N_{o}N_{s})^{2}$ . In terms of spin and orbital indices, the basis for the tensor on this form is given as all combinations of two indices $(\sigma\sigma^{\prime})\times(\alpha_{s}\beta_{s^{\prime}})$ . The spin combinations are $(\sigma\sigma^{\prime})=\{(\uparrow\uparrow),(\uparrow\downarrow),(\downarrow\uparrow),(\downarrow\downarrow)\}$ , and $(\alpha_{s}\beta_{s^{\prime}})$ are the 36 combinations of orbitals $\alpha=yz,xz,xy$ and sublattices $s=A,B$ .

Interactions are treated by calculating the irreducible particle-particle vertex, which is modified by spin-fluctuations from the random phase approximation (RPA). The irreducible vertex is given by the parquet equations [Bickers2004], where here the bare vertex is $\hat{\Gamma}_{0}=\hat{V}$ . The bare vertex is the Hubbard-Kanamori interaction [Kanamori1963], which is defined in real space as

$\displaystyle H_{\textrm{I}}=$	$\displaystyle U\sum_{\boldsymbol{j},\alpha}n_{\boldsymbol{j}\alpha\uparrow}n_{\boldsymbol{j}\alpha\downarrow}$	(18)
	$\displaystyle+\sum_{\boldsymbol{j},\alpha\neq\beta}J_{\textrm{H}}\left[c^{\dagger}_{\boldsymbol{j}\alpha\uparrow}c^{\dagger}_{\boldsymbol{j}\beta\downarrow}c_{\boldsymbol{j}\alpha\downarrow}c_{\boldsymbol{j}\beta\uparrow}+c^{\dagger}_{\boldsymbol{j}\alpha\uparrow}c^{\dagger}_{\boldsymbol{j}\alpha\downarrow}c_{\boldsymbol{j}\beta\downarrow}c_{\boldsymbol{j}\beta\uparrow}\right]$
	$\displaystyle+\sum_{\boldsymbol{j},\alpha<\beta,\sigma}\left[U^{\prime}n_{\boldsymbol{j}\alpha\sigma}n_{\boldsymbol{j}\beta\bar{\sigma}}+\left(U^{\prime}-J_{\textrm{H}}\right)n_{\boldsymbol{j}\alpha\sigma}n_{\boldsymbol{j}\beta\sigma}\right]$

with the intraorbital Hubbard interaction $U$ , the Hund’s coupling $J_{\textrm{H}}$ , and the interorbital repulsion $U^{\prime}=U-2J_{\textrm{H}}$ . Following the notation of Ref. [14], we define the bare vertex $\hat{V}$ as:

H_{\textrm{I}}=\displaystyle\frac{1}{N^{2}_{k}}\sum_{\boldsymbol{k},\boldsymbol{q}}\sum_{a,b,c,d}V_{abcd}(\boldsymbol{q})c^{\dagger}_{\boldsymbol{k},a}c^{\dagger}_{-\boldsymbol{k},c}c_{-(\boldsymbol{k}-\boldsymbol{q}),b}c_{\boldsymbol{k}-\boldsymbol{q},d}

(19)

for all orbit-spin-sublattice indices $a=(\alpha,\sigma,s)$ . On the same rank 4 tensor structure form as the susceptibility, the tensor $\hat{V}$ can be divided into spin sectors

\hat{V}=\left(\begin{array}[]{c c c c}\hat{V}^{\uparrow\uparrow\uparrow\uparrow}&0&0&\hat{V}^{\uparrow\uparrow\downarrow\downarrow}\\ 0&\hat{V}^{\uparrow\downarrow\uparrow\downarrow}&0&0\\ 0&0&\hat{V}^{\downarrow\uparrow\downarrow\uparrow}&0\\ \hat{V}^{\downarrow\downarrow\uparrow\uparrow}&0&0&\hat{V}^{\downarrow\downarrow\downarrow\downarrow}\end{array}\right)

(20)

with only on-site terms $s=s^{\prime}$ . Each matrix $\hat{V}^{\sigma_{1}\sigma_{2}\sigma_{3}\sigma_{4}}$ in the orbit-sublattice basis is defined as

V^{\uparrow\uparrow\uparrow\uparrow}_{\alpha_{s}\beta_{s}\gamma_{s}\delta_{s}}=V^{\downarrow\downarrow\downarrow\downarrow}_{\alpha_{s}\beta_{s}\gamma_{s}\delta_{s}}=\begin{cases}U^{\prime}-J_{\textrm{H}}&(\alpha=\gamma\neq\beta=\delta)\\ -U^{\prime}+J_{\textrm{H}}&(\alpha=\beta\neq\gamma=\delta)\\ 0&(\textrm{otherwise})\end{cases}

(21)

V^{\uparrow\uparrow\downarrow\downarrow}_{\alpha_{s}\beta_{s}\gamma_{s}\delta_{s}}=V^{\downarrow\downarrow\uparrow\uparrow}_{\alpha_{s}\beta_{s}\gamma_{s}\delta_{s}}=\begin{cases}-U&(\alpha=\beta=\gamma=\delta)\\ -J_{\textrm{H}}&(\alpha=\gamma\neq\beta=\delta)\\ -U^{\prime}&(\alpha=\beta\neq\gamma=\delta)\\ -J_{\textrm{H}}&(\alpha=\delta\neq\beta=\gamma)\\ 0&(\textrm{otherwise})\end{cases}

(22)

V^{\uparrow\downarrow\uparrow\downarrow}_{\alpha_{s}\beta_{s}\gamma_{s}\delta_{s}}=V^{\downarrow\uparrow\downarrow\uparrow}_{\alpha_{s}\beta_{s}\gamma_{s}\delta_{s}}=\begin{cases}U&(\alpha=\beta=\gamma=\delta)\\ U^{\prime}&(\alpha=\gamma\neq\beta=\delta)\\ J_{\textrm{H}}&(\alpha=\beta\neq\gamma=\delta)\\ J_{\textrm{H}}&(\alpha=\delta\neq\beta=\gamma)\\ 0&(\textrm{otherwise})\end{cases}

(23)

for each sublattice $s=A,B$ and for orbital indices $\alpha,\beta,\gamma,\delta$ . The effective particle-particle vertex [Bickers2004] is:

\hat{\Gamma}^{pp}(q=K-K^{\prime})=-\frac{1}{2}\hat{V}-\hat{V}\hat{\chi}(K-K^{\prime})\hat{V}

(24)

The susceptibility $\hat{\chi}(q)$ is approximated as the RPA susceptibility:

\hat{\chi}(q)=\left(1-\hat{\chi}_{0}(q)\hat{V}\right)^{-1}\hat{\chi}_{0}(q)

(25)

which uses the same bare vertex $\hat{V}$ and has the same tensor structure as $\hat{\chi}_{0}(q)$ . The tensor multiplication is defined as a matrix multiplication in the spin and orbital-sublattice combinations

\left[\hat{A}\hat{B}\right]^{(\sigma\sigma^{\prime})_{1}(\sigma\sigma^{\prime})_{2}}_{(\alpha_{s}\beta_{s^{\prime}})_{1}(\alpha_{s}\beta_{s^{\prime}})_{2}}=\sum_{\tilde{\alpha},\tilde{\beta},\tilde{\sigma},\tilde{\sigma}^{\prime}}A^{(\sigma\sigma^{\prime})_{1}(\tilde{\sigma}\tilde{\sigma}^{\prime})}_{(\alpha_{s}\beta_{s^{\prime}})_{1}(\tilde{\alpha}\tilde{\beta})}B^{(\tilde{\sigma}\tilde{\sigma}^{\prime})(\sigma\sigma^{\prime})_{2}}_{(\tilde{\alpha}\tilde{\beta})(\alpha_{s}\beta_{s^{\prime}})_{2}}

(26)

The effective vertex is not further spin-diagonalized into spin-singlet and spin-triplet vertices as our regime of intermediate to strong spin-orbit coupling inherently mixes the two sectors. The linearized gap equation can be obtained from the Luttinger-Ward functional [Baym1961, Baym1962], where the linearization of the anomalous Green’s function is given by the Dyson-Gorkov equations Eq. (9).

\Delta_{ab}(k)=\frac{1}{N\beta}\sum_{k^{\prime}}\sum_{a^{\prime}b^{\prime}}\Gamma^{pp}_{aa^{\prime}b^{\prime}b}(k-k^{\prime})F_{a^{\prime}b^{\prime}}(k^{\prime})

(27)

F_{a^{\prime}b^{\prime}}(k^{\prime})=\sum_{\mu\nu}G_{a^{\prime}\mu}(k^{\prime})\bar{G}_{\nu b^{\prime}}(k^{\prime})\Delta_{\mu\nu}(k^{\prime})

(28)

where $\beta$ is the inverse temperature and $N=N_{k}^{2}$ . Each index here runs over all orbit-spin-sublattice combinations $a=(\alpha,\sigma,s)$ . Further, we apply both the static ( $\omega_{n}\rightarrow\delta=10^{-4}$ eV, $\hat{\Delta}(k)\rightarrow\hat{\Delta}(\boldsymbol{k})$ , $\Gamma^{pp}(q)\rightarrow\Gamma^{pp}(\boldsymbol{q})$ ) and normal state ( $\hat{G}(k^{\prime})\rightarrow\hat{G}^{(0)}(k^{\prime})$ ) approximations. The linearized gap equation can be solved as an eigenvalue problem, as a version of the Eliashberg equation:

\lambda_{e}\Delta_{ab}(\boldsymbol{k})=\frac{1}{N}\sum_{\boldsymbol{k}^{\prime},a^{\prime}b^{\prime},\mu\nu}\Gamma^{pp}_{aa^{\prime}b^{\prime}b}(\boldsymbol{k}-\boldsymbol{k}^{\prime})\phi_{a^{\prime}b^{\prime}\mu\nu}^{\boldsymbol{k}^{\prime}}\Delta_{\mu\nu}(\boldsymbol{k}^{\prime})

(29)

with

\phi_{a^{\prime}b^{\prime}\mu\nu}^{\boldsymbol{k}^{\prime}}=-\frac{1}{\beta}\sum_{\omega_{n}}G^{(0)}_{a^{\prime}\mu}(k^{\prime})\bar{G}^{(0)}_{\nu b^{\prime}}(k^{\prime})

(30)

A largest eigenvalue of unity or higher $\lambda_{e}\geq 1$ indicates a possible superconducting order. A non-explicit summation over Matsubara frequencies is performed, like in Eq. (17), and the equation depends only on momentum $\boldsymbol{k}$ .

While this system is strongly interacting one might expect a multitude of order parameters including magnetic orders. We use the Stoner criterion to identify phases in the particle-hole channel. Methods attempting to treat particle-hole and particle-particle self-energies on equal footing are beyond the scope of this work [Bickers1989, Bickers1991]. The RPA susceptibility, Eq. (25), will pass a critical point and diverge, when the Hartree-Fock term in the particle-hole channel has an eigenvalue above unity at any $k$ -point which we name $Q$ [Bickers1989, Bickers1991]:

\textrm{max eig}\left[\hat{\chi}_{0}(Q)\hat{V}\right]=1

(31)

The particle-hole instability is given by the eigenvector to the tensor $\hat{\chi}_{0}(Q)\hat{V}$ , which can be unfolded into a matrix. A rank 4 tensor $C_{ijkl}$ of dimension $N_{1}\times N_{2}\times N_{1}\times N_{2}$ can be mapped onto a matrix $C_{\mu\nu}$ of dimension $N_{1}N_{2}\times N_{1}N_{2}$ via

$\displaystyle\mu=i+N_{1}(j-1),$	$\displaystyle\qquad\mu\in\left[1,\dots,N_{1}N_{2}\right]$	(32)
$\displaystyle i=1+\textrm{mod}(\mu-1,N_{1}),$	$\displaystyle\qquad i\in\left[1,\dots,N_{1}\right]$
$\displaystyle j=1+\textrm{div}(\mu-1,N_{1}),$	$\displaystyle\qquad j\in\left[1,\dots,N_{2}\right]$

where ”mod” is the modulus operation and ”div” is integer division. The given mapping preserves the defined tensor multiplication, in Eq. (26), as matrix multiplication in the unfolded matrix. The type of instability can also be classified by the magnetic channel it occurs in. The susceptibility can be spin block-diagonalized by rewriting it in the basis of magnetic operators $m^{z}_{\tilde{\alpha}\tilde{\beta}}=\frac{1}{\sqrt{2}}(c^{\dagger}_{\tilde{\alpha}\uparrow}c_{\tilde{\beta}\uparrow}-c^{\dagger}_{\tilde{\alpha}\downarrow}c_{\tilde{\beta}\downarrow})$ :

\hat{\chi}^{z}_{\tilde{\alpha}\tilde{\beta}\tilde{\gamma}\tilde{\delta}}(q)=\frac{1}{N_{k}^{2}}\int_{0}^{\beta}d\tau e^{i\omega_{n}\tau}\langle T_{\tau}m^{z}_{\tilde{\alpha}\tilde{\beta},\boldsymbol{q}}(\tau)m^{z}_{\tilde{\gamma}\tilde{\delta},-\boldsymbol{q}}(0)\rangle_{c}

(33)

with the indices $\tilde{\alpha}=\alpha_{s}$ running over all orbital $\alpha=yz,xz,xy$ and sublattice $s=A,B$ combinations. The magnetic channels are

\hat{\chi}^{z}(q)=\frac{1}{2}\left(\hat{\chi}^{\uparrow\uparrow\uparrow\uparrow}(q)-\hat{\chi}^{\uparrow\uparrow\downarrow\downarrow}(q)\right)

(34)

\hat{\chi}^{+}(q)=\hat{\chi}^{\uparrow\downarrow\uparrow\downarrow}(q),\qquad\hat{\chi}^{-}(q)=\hat{\chi}^{\downarrow\uparrow\downarrow\uparrow}(q)

(35)

with the out-of-plane spin $\hat{\chi}^{z}(q)$ and in-plane spin $\hat{\chi}^{\pm}(q)$ channels ¹¹1The spin block-diagonalization introduces an additional density channel $\hat{\chi}^{d}(q)$ . No density-channel instability is found in this work and the peaks are significantly smaller than in the magnetic channels. This remains true for $\hat{\chi}^{d}_{J}(q)$ in the $j$ -state basis.. As the $j$ -state basis describes the bands better than the orbital-spin basis, the $\chi(q)$ -terms can be transformed, via Eq. (8), as

\chi_{J,ijkl}(q)=\sum_{\alpha\beta\gamma\delta}M_{i\alpha}M^{\dagger}_{\gamma k}M_{l\delta}M^{\dagger}_{\beta j}\chi_{\alpha\beta\gamma\delta}(q)

(36)

where $M$ is the transformation from spin and orbit to the total angular momentum basis. In this basis, the susceptibility is similarly divided into different pseudospin channels $\hat{\chi}^{z}_{J}(q)$ and $\hat{\chi}^{\pm}_{J}(q)$ .

II.3 Computational details

All RPA calculations, for the susceptibility Eq. (25) and superconductivity Eq. (29) were performed on a $N_{k}\times N_{k}=46\times 46$ momentum lattice. The finite momentum resolution limits the lowest accessible temperature. Peaks in the susceptibility cannot be narrower than the lattice spacing and we thus require some thermal broadening to get reliable results. In this work we use $k_{\textrm{B}}T=0.001$ eV ( $T\approx 11$ K). The temperature was chosen such that less than a $10\%$ change in largest eigenvalue was found when going from a lattice of size $N_{k}=32$ to $N_{k}=46$ , at most points. In the data of Figs. 7b & 11 a larger change is observed, therefore calculations for these figures where done at $N_{k}=64$ . The specific temperature is of interest for a potential $d$ -wave superconducting order in the electron doped compound. The 2016 experiment in Ref. [16] observed an order with this symmetry below $T\approx 30$ K, with a maximum at $T\approx 10$ K. The largest value of the linearized gap equation is found via the Arnoldi method, with a convergence criterion of $10^{-7}$ .

III Phase diagram: strain and doping

In the following sections we discuss the phase diagrams obtained by varying the charge doping and applying an increasing compressive strain. For the most realistic parameter range one obtains Fig. 3. To explore additional effects from chemical doping and to illustrate the richness of similar compounds, additional phase diagrams are shown in Fig. 4. All phase diagrams show strain-induced superconductivity for a broad range of parameters. In addition, three main features can be observed. First, the two types of superconductivity found in earlier works [11, 14], a pseudospin $j=1/2$ $d$ -wave and an orbital $s_{\pm}$ -wave, are found. In addition, both types can be found when varying only the doping, for a select set of parameters. Second, at high enough strain an orbital selective pairing which favors one of the in-plane directions can be found. This type of order has a larger component originating from the $xz$ orbital than from $yz$ , spontaneously breaking the in-plane symmetry. And third, under compressive strain the susceptibility goes from having the largest peaks in the pseudospin states to peaks instead originating from the spin states. This shift affects both the magnetic order and the possible pairing symmetries.

Choosing a realistic regime for the phase diagram, there are two criteria for the chosen interaction parameters. The first criterion is that sufficient doping, either hole or electron, should in accordance to experiment, result in a transition out of the magnetic order. The second criterion is that in undoped Sr₂IrO₄ the magnetic order persists up to a strain value of $\epsilon\approx-2\%$ [31, 39]. For the most realistic values of the Hund’s and spin-orbit coupling, $J_{\textrm{H}}=0.1U$ and $\lambda=0.6$ eV, the first criterion is satisfied for $U\approx 2|t|$ , shown in Appendix A. The realistic phase diagram is presented for $U=1.1$ eV $\approx 3|t|$ . Since the model overestimates the orders this choice only satisfies the second criterion. In Appendix A the phase diagram for $U=1.4$ eV $\approx 4|t|$ results in the same phase transitions at higher strain values. Complementary mean field calculations, containing magnetic, superconducting as well as other order parameters, yields a qualitatively similar phase diagram, with differences explained in Appendix B.

In Fig. 4, we also consider a higher Hund’s coupling of $J_{\textrm{H}}=0.25U$ , with a lower spin orbit coupling of $\lambda=0.5$ eV for the following reason. Hole doping via the substitution of iridium atoms for rhodium or ruthenium atoms could modify the effective interaction parameters and spin-orbit coupling. Some works estimate the spin orbit coupling of iridium, rhodium and ruthenium as $\lambda^{\textrm{Ir}}\approx 0.45$ eV, $\lambda^{\textrm{Rh}}=\lambda^{\textrm{Ru}}\approx 0.19$ eV respectively [Lee2021, Zwartsenberg2020, Brouet2021]. Moreover, ruthenium atoms have a higher Hund’s coupling of $J_{\textrm{H}}\approx 0.15U-0.2U$ [Georges2013]. The full set of phase diagrams are presented in Figs. 3 & 4, with identified superconducting and magnetic phases.

IV Magnetic order

In the RPA calculation, the particle-hole order is found using the Stoner criterion, Eq. (31). An order can be characterized by two features of $\hat{\chi}(\boldsymbol{q})$ . First is the nesting vector $\boldsymbol{Q}=(q_{x},q_{y})$ , given in the extended BZ, at which the Stoner criterion is met. The instabilities in the phase diagrams all appear at the four points, $\boldsymbol{Q}_{1}\approx(\pi,\pi)$ , $\boldsymbol{Q}_{2}\approx(\frac{\pi}{2},\frac{\pi}{2})$ , $\boldsymbol{Q}_{0}\approx(0,0)$ , and $\boldsymbol{Q}_{1d}\approx(\frac{\pi}{2},\pi)$ . The exact location of the instabilities is shifted a small distance from the ideal values, which depends on doping as well as strain. Examples of the dominant peaks are shown in Fig. 5. The most relevant components of the susceptibility are shown in Figs. 6 & 7, along a few cuts in the phase diagrams. For the lower Hund’s coupling in Fig. 3 the only instability is at $\boldsymbol{Q}_{1}$ . An order described in real space by a two site unit cell in a square lattice will have a reduced Brillouin zone with a unit vector $\boldsymbol{Q}=(\pi,\pi)$ . As antiferromagnetism is a two site order expected in this compound, the $\boldsymbol{Q}_{1}$ susceptibility peaks are expected to be antiferromagnetic. In Fig. 4 the $\boldsymbol{Q}_{2}$ -instability occurs at most doping values. A real space order corresponding to $\boldsymbol{Q}=(\frac{\pi}{2},\frac{\pi}{2})$ , will have a unit cell containing 4 sites. However, it should be noted that the peak $\boldsymbol{Q}_{2}\approx(\frac{\pi}{2},\frac{\pi}{2})$ is doping dependent and never occurs exactly at this value. It is an incommensurate order closer to $(0.6\pi,0.6\pi)$ for hole doping and $(0.4\pi,0.4\pi)$ for electron doping. In Fig. 5 there is an additional copy of these peaks that connects different copies of pockets in the extended BZ. Another instability, at $\boldsymbol{Q}_{1d}$ , only occurs at the highest strains and electron dopings considered. This nesting vector connects segments on the FS with either a clear $xz$ -character to other segments belonging to the same orbital. The quasi-1d dispersion of the $xz$ orbital, with $t_{1}\gg t_{\delta}$ in Eq. (4), leads to the susceptibility in the intra- $xz$ channel having peaks connecting point in momentum space mainly along the $x$ -direction. Finally, a ferromagnetic instability at $\boldsymbol{Q}_{0}$ is possible in Fig. 4a at $n=4.7$ .

The second feature that characterizes a magnetic instability is the channel in which the instability occurs. In Eqs. 34 & 35, the spin susceptibility $\hat{\chi}(\boldsymbol{q})$ , as well as the pseudospin susceptibility $\hat{\chi}_{J}(\boldsymbol{q})$ , are divided into magnetic channels. Most instabilities found in Figs. 3 & 4 are of the type $\chi^{\uparrow\downarrow\uparrow\downarrow}_{J,\left(\frac{1}{2},\pm\frac{1}{2}\right)}(\boldsymbol{Q}_{1})$ , an in-plane magnetic instability with mainly $j=1/2$ contributions. This can be interpreted as the canted in-plane antiferromagnetic order (x-cAFM) observed in each layer experimentally. As denoted in Fig. 6b, the higher Hund’s coupling, $J_{\textrm{H}}=0.25U$ and $U=1.1$ eV, results in an additional out-of-plane component accompanying the in-plane order, with $\chi^{z}_{J,\left(\frac{1}{2},\pm\frac{1}{2}\right)}(\boldsymbol{Q}_{1})$ .

At any point where the $\boldsymbol{Q}_{2}$ instability is present, it occurs in channels of the spin susceptibility, rather than pseudospin. The instability is in-plane and has contributions mainly from the $yz$ and $xz$ orbitals: $\chi^{\uparrow\downarrow\uparrow\downarrow}_{yz}(\boldsymbol{Q}_{2})$ & $\chi^{\uparrow\downarrow\uparrow\downarrow}_{xz}(\boldsymbol{Q}_{2})$ . Even though this instability is only present at $U=1.4$ eV and $J_{\textrm{H}}=0.25U$ , these susceptibility peaks remain large in the entire Fig. 4a phase diagram. At high hole doping in Fig. 4a a purely out-of-plane ferromagnetic instability $\chi^{z}(\boldsymbol{Q}_{0})$ occurs for spins in each of the orbitals $yz$ , $xz$ , and $xy$ .

As the strain increases the peaks in the susceptibility at different channels decrease at different rates. Moreover, we will see below that superconductivity found directly adjacent to a magnetic order is mediated by those magnetic fluctuations while superconducting orders which are mediated by other fluctuations can become more favorable as strain is increased further. In Fig. 6b, at $J_{\textrm{H}}=0.25U$ , the spin susceptibility peaks decrease significantly slower than those of the antiferromagnetic pseudospin order. At strain $\epsilon=-0.5\%$ , and beyond, the peaks in spin susceptibility become dominant. The effective particle-particle vertex, Eq. (24), develops peaks at the same $\boldsymbol{Q}$ -points as for the spin (or pseudospin) susceptibility. To track which magnetic fluctuation mediates the superconducting orders, the strengths of spin and pseudospin fluctuations are compared. The total susceptibilities in the two bases are

\chi^{\sigma\sigma^{\prime}\sigma\sigma^{\prime}}(\boldsymbol{q})=\frac{1}{2}\sum_{s,s^{\prime}}\sum_{\alpha,\beta}e^{iq_{x}(\Theta_{s}-\Theta_{s^{\prime}})}\chi^{(\alpha\sigma)(\alpha\sigma^{\prime})(\beta\sigma)(\beta\sigma^{\prime})}_{sss^{\prime}s^{\prime}}(\boldsymbol{q})

(37)

\chi^{\tau\tau^{\prime}\tau\tau^{\prime}}_{J}(\boldsymbol{q})=\frac{1}{2}\sum_{s,s^{\prime}}\sum_{m,n}e^{iq_{x}(\Theta_{s}-\Theta_{s^{\prime}})}\chi^{(m\tau)(m\tau^{\prime})(n\tau)(n\tau^{\prime})}_{J,sss^{\prime}s^{\prime}}(\boldsymbol{q})

(38)

where $\Theta_{s}$ , with $s=A,B$ , chooses the sublattice such that $\Theta_{A}=0$ and $\Theta_{B}=1$ .

It should be noted, that the pseudospin susceptibility nesting vector $\boldsymbol{Q}_{1}$ can be related to the hidden spin density wave (hSDW) orders found in other models for multi-band superconductors [Berg2012, Rodriguez2020]. The $j$ -state basis mixes spin and orbital degrees of freedom and therefore the peaks found correspond to a linear combination of channels that favors both SDW and hSDW orders.

V Superconductivity

Symmetry $\eta^{\mu}_{R}(\boldsymbol{k})$	$R=0$	$R=1$	$R=2$	$R=3$
A_1g, $s$	1	$\cos k_{x}+\cos k_{y}$	$2\cos k_{x}\cos k_{y}$	$\cos 2k_{x}+\cos 2k_{y}$
B_1g, $d_{x^{2}-y^{2}}$	-	$\cos k_{x}-\cos k_{y}$	-	$\cos 2k_{x}-\cos 2k_{y}$
B_2g, $d_{xy}$	-	-	$2\sin k_{x}\sin k_{y}$	-
E_u, $p$	-	$\sin k_{x},\sin k_{y}$	$\sin(k_{x}+k_{y}),$	$\sin 2k_{x},\sin 2k_{y}$
			$\sin(k_{x}-k_{y})$

Table 1: The lattice harmonics for each relevant symmetry on the square lattice is given for different radii

R

, describing what distance neighbors the symmetry is found on.

R=0

is on-site,

R=1

is nearest neighbors, and so on.

V.1 Symmetries

In single-orbital models, it is of highest importance to determine whether the superconductivity is a spin-singlet or a spin-triplet order. Topological superconductivity and Majorana modes arise from superconductivity with $p$ -wave pairing, which requires spin-triplet pairing in those systems. Efforts to induce superconductivity through the proximity effect are thus often focused on finding spin-triplet orders. However, once multiple orbitals and spin-orbit coupling are considered the connection between spin-triplets and $p$ -wave symmetry is no longer a strict requirement. Multi-orbital models allow for a large set of possible pairing symmetries. The pairing matrix $\hat{\Delta}(\boldsymbol{k})$ must be antisymmetric under the full $\mathcal{S}\mathcal{P}\mathcal{O}\mathcal{T}$ -exchange [Sigrist1991, Linder2019]. Therefore, any pairing can be classified as being either even or odd under spin exchange ( $\mathcal{S}$ ), relative coordinate reflection ( $\mathcal{P}$ ), orbital exchange ( $\mathcal{O}$ ), and relative time exchange ( $\mathcal{T}$ ) as defined in Appendix D. As only the static case is considered in this work the order parameter is constant, and therefore even, under $\mathcal{T}$ . Note that the operators $\mathcal{P}$ and $\mathcal{T}$ only exchange relative parameters, and are thus different from the reflection and time reversal operators.

For example, classifying symmetries for strong spin-orbit coupling in the predicted $j=1/2$ $d$ -wave, will inherently result in both spin-singlet and spin-triplet contributions. The pairing is more accurately described by the total angular momentum of the pair, which can take the values $J=0,1,2,3$ [Dutta2021, Adarsh2022]. Only within the $j_{1}\otimes j_{2}=\frac{1}{2}\otimes\frac{1}{2}$ sector do we still only get pairs that are either a $J=0$ singlet or triplets $J=1$ with $M=-1,0,+1$ . In Appendix D the symmetry operations in the orbital basis are shown for the two types of pairing found in Fig. 3.

The spatial symmetry can be considered for the non-interacting bands $n_{\tau}$ by projecting the intraband pairing onto the Fermi surface (FS). The FS may contain three types of pockets belonging to two types of bands. The larger pockets, consisting mainly of $(j,j_{z})=(\frac{1}{2},\pm\frac{1}{2})$ states, are located around the points $(k_{1},k_{2})=(\pi,\pi),(\pi,0)$ and have superconducting gaps that can be described by the same spatial symmetry. We therefore only look at one of these pockets, denoted FS₁. The smaller FS₂, with mainly $(j,j_{z})=(\frac{3}{2},\pm\frac{3}{2})$ is centered around $(k_{1},k_{2})=(0,0)$ . The spatial symmetry is thus studied for four intraband parameters: pseudospin-singlets $\Delta^{s}_{\textrm{FS}_{1}}(\boldsymbol{k})$ , $\Delta^{s}_{\textrm{FS}_{2}}(\boldsymbol{k})$ and pseudospin-triplets $\Delta^{t}_{\textrm{FS}_{1}}(\boldsymbol{k})$ , $\Delta^{t}_{\textrm{FS}_{2}}(\boldsymbol{k})$ . $\Delta^{s/t}_{\textrm{FS}_{n}}(\boldsymbol{k})$ is the order parameter Eq. (29) projected onto the band at the Fermi surface $FS_{n}$ . The number of points belonging to a pocket is $N_{\textrm{FS}_{n}}$ . For example the pseudospin-singlet is

$\displaystyle\Delta^{s}_{\textrm{FS}_{n}}(\boldsymbol{k})=$	$\displaystyle\sum_{\boldsymbol{k}\in\textrm{FS}_{n}}\sum_{\alpha\beta}\sum_{mm^{\prime}}\left[U_{\boldsymbol{k}}^{\dagger}\right]_{n\alpha}\left[U_{-\boldsymbol{k}}^{\dagger}\right]_{n\beta}$	(39)
	$\displaystyle\times\frac{1}{\sqrt{2}}\left(M_{\alpha(m+)}^{\dagger}M_{\beta(m^{\prime}-)}^{\dagger}\Delta_{(m+)(m^{\prime}-)}(\boldsymbol{k})\right.$
	$\displaystyle\left.-M_{\alpha(m-)}^{\dagger}M_{\beta(m^{\prime}+)}^{\dagger}\Delta_{(m-)(m^{\prime}+)}(\boldsymbol{k})\right)$

where $M$ is the matrix in Eq. (8) and the matrix $U_{\boldsymbol{k},\alpha n}$ transforms the band basis into the orbital basis. In the numerical calculations, there are small but non-zero interband contributions that will not be considered further. The spatial symmetries can be quantified via the projection coefficients $P^{\mu,R}_{s/t,\textrm{FS}_{n}}$ onto the basis functions for each parity irreducible representation, $\eta^{\mu}_{R}(\boldsymbol{k})$ , as given in Table 1

\Delta^{s/t}_{\textrm{FS}_{n}}(\boldsymbol{k})=\sum_{\mu,R}P^{\mu,R}_{s/t,\textrm{FS}_{n}}\eta^{\mu}_{R}(\boldsymbol{k}).

(40)

The projection coefficient for each irreducible representation is thus found via

P^{\mu,R}_{s/t,\textrm{FS}_{n}}=\frac{1}{N_{\textrm{FS}_{n}}}\sum_{\boldsymbol{k}\in\textrm{FS}_{n}}\eta^{\mu}_{R}(\boldsymbol{k})\Delta^{s/t}_{\textrm{FS}_{n}}(\boldsymbol{k}).

(41)

Examples of the two leading types of pairings are shown in Fig. 8. The most prominent $d$ -wave is a $j=1/2$ pseudospin singlet with a $\eta^{B_{1g}}_{R=1}$ pairing. The found $s_{\pm}$ -wave order has contributions from both $\eta^{A_{1g}}_{R=0}(\boldsymbol{k})$ and $\eta^{A_{1g}}_{R=2}(\boldsymbol{k})$ , with opposite sign for the two pockets. The main $j$ -states components are from $(\frac{1}{2},\pm\frac{1}{2})$ and $(\frac{3}{2},\pm\frac{3}{2})$ . However, components mixing $j$ -states is stronger than for the $d$ -wave. The symmetry in terms of orbital origin is specified further in Appendix D.

V.2 Realistic strain-induced order

At low Hund’s coupling and high SOC, a $d$ -wave is found for a wide range of doping values once there is no longer a magnetic order present. The mean field calculations in Appendix B corroborate the prediction of this order. The magnetic region extends up to $\epsilon=-1.5\%$ at the undoped $n=5$ , while it persists at higher strains on the electron doped side. Because of the required compressive strain, two bands are present at the FS at all points where superconductivity is found. As seen in Fig. 8, the gap on FS₁ is significantly larger than on FS₂. The $d$ -wave originates from the $j=1/2$ states, which are the states the band at FS₁ also belongs to. Strain has increased the size of the $j=1/2$ electron pocket, FS₁, in the entire region where superconductivity is found. Even for hole doping at $\epsilon=-2\%$ the pocket FS₁ is comparable in size to FS₁ at $\epsilon=0$ and electron doping, as further explained in Appendix B.

V.3 Varying Hund’s coupling

There are several important factors that determine which pairing symmetry is favored. Doping and strain change the shape of the Fermi surface, and thus the nesting vectors $\boldsymbol{Q}$ , as well as the orbital contributions in each pocket. Both compressive strain and hole doping increase the presence of $yz$ and $xz$ orbitals. However, the type of fluctuations which dominate the RPA interaction vertex depends on the interaction parameters. In Fig. 8 the largest eigenvalue of Eq. (29), and the symmetry of the pairing, are shown as the Hund’s coupling varies. We observe a general trend in which the $s_{\pm}$ -order becomes more favorable than $d$ -wave superconductivity at $J_{\textrm{H}}\geq 0.25U$ . For the lower SOC, $\lambda=0.5$ eV, the largest eigenvalue is above unity for all values of Hund’s coupling considered. The $s_{\pm}$ -wave is only possible for low SOC and hole doping, since that places the Fermi level deeper in the band of $(j,j_{z})=(\frac{3}{2},\pm\frac{3}{2})$ character. By contrast, the $d$ -wave order is present for a wider range of parameters.

At all points there are two main types of competing fluctuations; pseudospin fluctuations around $\boldsymbol{Q}_{1}\approx(\pi,\pi)$ and spin fluctuations around $\boldsymbol{Q}_{2}\approx(\frac{\pi}{2},\frac{\pi}{2})$ . By studying the maximum peak values in Fig. 9 one can determine which one of these fluctuations best describes the system. As can be seen in Figs. 6 & 7, dominating spin fluctuations promote an $s_{\pm}$ -wave order. However, if both types of fluctuations are of roughly equal size the $j=1/2$ $d$ -wave is favored.

V.4 Multi-band orders

For a higher Hund’s coupling superconducting orders are favored which open large gaps on multiple pockets at the Fermi surface. In Fig. 4a, there are two distinct superconducting orders at $U=1.1$ eV. The magnetic order disappears for small strains and superconductivity is only possible up to $\epsilon=-2\%$ . Moreover, the symmetry of the superconducting order is dependent on doping. At high hole doping and some strain the pairing is an $s_{\pm}$ -wave, while remaining a $d$ -wave for all other doping values. When we set the Hubbard $U$ to a higher value of $U\approx 4|t|$ , as seen in Fig. 4, and turn on a high compressive strain the band structure changes. This can lead to new pairing functions which were not seen earlier. A drastic change occurs in Fig. 4b where some regions have an anisotropic $s$ - or $d$ -wave pairing. For all values of $J_{\textrm{H}}$ considered in Fig. 10, the anisotropic order, which is a mix of $s$ - and $d$ - wave pairing, is found. The $s$ -wave contribution increases with Hund’s. The orbital components, as well as a simple model for this state, are described in Appendix E. This superconducting order is an orbital-selective state with a stronger spin-singlet in the $xz$ -orbital. One notable reason for this new type of pairing is the increased spin nature of the fluctuations as compressive strain is increased. As already discussed, all regions with anisotropic pairing are mediated by a large spin susceptibility peak in $\chi^{\uparrow\downarrow\uparrow\downarrow}_{xz}(\boldsymbol{Q}_{1d})$ . In Fig. 9, the larger spin to pseudospin susceptibility ratio can be compared for the increased strain, at all $J_{\textrm{H}}$ values.

V.5 Odd parity

In general, the out-of-plane $\chi^{z}(\boldsymbol{Q}_{0})$ susceptibility peaks remain smaller than either the in-plane spin $\chi^{\uparrow\downarrow\uparrow\downarrow}\left(\boldsymbol{Q}\right)$ or pseudospin $\chi^{\uparrow\downarrow\uparrow\downarrow}_{J}\left(\boldsymbol{Q}\right)$ peaks. This is the case for all values calculated so far, except for the one small region in Fig. 7 with large hole doping, high Hund’s coupling, and low compressive strain. A high enough Hubbard interaction $U\geq 1.1$ eV is also required. As several of these factors also increase the in-plane spin susceptibility, the out-of-plane susceptibility only dominates in a very small parameter range. In Fig. 11 one such small patch can be found for a very high Hund’s coupling, $J_{\textrm{H}}=0.3U$ . Accompanying these fluctuations is an odd parity $p$ -wave superconductivity. The main contribution, described by the symmetry representation detailed in Appendix D, is

\Delta_{p}(\boldsymbol{k})\approx\left(h_{x}\otimes\sigma^{z}+ih_{y}\otimes\mathbb{I}\right)\otimes\lambda_{3}.

(42)

where $\lambda_{j}$ is the $j^{th}$ Gell-Mann matrix [GellMann1961] acting in the three-orbitals space. The $p$ -wave pairing has several leading terms, which are $h_{x}=-\sin k_{x}+\sin k_{x}\cos k_{y}+\sin(2k_{x})$ and $h_{y}=-\sin k_{y}+\sin k_{y}\cos k_{x}+\sin(2k_{y})$ . Some smaller terms are proportional to $\left(h_{x}\otimes\sigma^{z}+ih_{y}\otimes\mathbb{I}\right)\otimes\lambda_{1}$ and $\Delta^{\downarrow\downarrow}_{xy,xy}(\boldsymbol{k})\propto\Delta^{\uparrow\uparrow}_{yz,yz}(\boldsymbol{k})$ . Up- and down-spin pairing have the opposite chirality. This helical $p$ -wave thus preserves time-reversal symmetry (TRS), since $\Delta_{\alpha\beta}^{\downarrow\downarrow}(\boldsymbol{k})=\left(\Delta_{\alpha\beta}^{\uparrow\uparrow}\right)^{\ast}(-\boldsymbol{k})$ . The helical nature is preserved for both Fermi surfaces, FS₁ and FS₂, with a larger gap on FS₂. A $\mathbb{Z}_{2}$ topological invariant can therefore be determined. We can consider a Chern number for each pseudospin sector, as outlined in Appendix F. The defining features of the pairing are

•

The helical $p$ -wave preserves time reversal symmetry, such that the total Chern number vanishes, $\mathcal{C}_{\textrm{tot}}=0$ .
•

The pseudospin Chern number, $\nu=(\mathcal{C}_{+}-\mathcal{C}_{-})/2$ , also vanishes such that in each pseudospin sector $\mathcal{C}_{\tau}=0$ .
•

In each pseudospin sector we define the pocket Chern number $\mathcal{C}_{n,\tau}$ for the band $n$ . We find $\mathcal{C}_{n,\tau}=\pm 1$ such that for each band the pseudospin Chern number is $\nu_{n}=\pm 1$ . The pocket contributions cancel such that $\nu=\sum_{n}\nu_{n}=0$ , as the pocket of $j=3/2$ - character has the opposite chirality to the $j=1/2$ pocket.
•

Therefore no topologically protected edge/vortex modes are expected.

Even for chiral TRS-breaking superconductors the Fermi surface topology can result in a topologically trivial state, when multiple pockets are present [Sato2009, Farrell2013]. The possible Chern numbers in multiband superconductors depend on the pairing function, and the location in the Brillouin zone of the resulting topological charges [Raghu2010, Imai2012, Scaffidi2015].

It might be possible to expand the $p$ -wave regime by tuning parameters such that ferromagnetism is favored. In our model, increasing the Hubbard coupling, $U$ , accomplishes this. However, with a higher $U$ , a higher strain required to reach the superconducting regime but a higher $U$ also increases the in-plane fluctuations, as in Fig. 7. This is caused by the other fluctuations being favored when the pocket with $j=3/2$ states becomes large. Only for a smaller pocket and a large $U$ would the $p$ -wave be favored. The odd parity order is thus a fine-tuned case which is found beyond realistic parameters.

VI Discussion

While strong interactions, spin-orbit coupling, and the proximity of multiple $d$ -bands to the Fermi level all point to the possibility of unconventional superconductivity, experimental observation of superconductivity is still missing. In this work we suggest that compressive strain may be a possible knob that, together with doping, can turn the system from magnetic to superconducting. We model the system using the extended Hubbard-Kanamori Hamiltonian and map out its phase diagram. Magnetic orders are found using the Stoner criterion while superconductivity is studied using the RPA linearized Eliashberg equation. For the range of parameters considered in this work, we find prominent regions of strain-induced superconductivity, among them a large fraction exhibits $d$ -wave pairing.

The $d$ - and $s_{\pm}$ -wave orders found are the same as in previous studies of the unstrained compound. For the values considered here the $s_{\pm}$ -wave can arise adjacent to the AFM order, for a high enough Hund’s coupling. For high strains and Hund’s coupling, new orbital-selective, anisotropic $s$ - or $d$ -wave orders are found. In addition, we find an out-of-plane ferromagnetic order. In a very fine-tuned region, the out-of-plane susceptibility mediates an odd parity $p$ -wave order. We note that the work of Ref. [12] finds a $p$ -wave order in hole doped Sr₂IrO₄, albeit at an extremely large Hubbard coupling, $U=12t$ , and a lower Hund’s coupling of $J_{\textrm{H}}=0.15U$ . The odd parity order we find is favored by a high Hubbard coupling, $U$ , in the unstrained compound and could be the dominant order at those values. At higher strains the ferromagnetic fluctuations vanish. The $p$ -wave order is found to be helical and topologically trivial, as determined via the $\mathbb{Z}_{2}$ invariant. However, other values of the Hund’s coupling and different Fermi surface geometry could potentially break time reversal symmetry and result in a chiral $p$ -wave [Scaffidi2014]. As in the case of Sr₂RuO₄, the topological nature of such a state in the iridates is highly dependent on the Fermi surface geometry and orbital composition.

It should be mentioned that the type of possible magnetic instabilities is not altered by the compressive strain for the realistic value of the Hund’s coupling, $J_{\textrm{H}}=0.1U$ . However, with increasing compression the pseudospin $j=1/2$ susceptibility decreases faster than the spin susceptibility, bringing the two leading fluctuation peaks to comparable sizes. At high Hund’s coupling and lower SOC we find an antiferromagnetic order that can have a mixed in- and out-of-plane structure. Moreover, at higher strains a spin-like incommensurate magnetic order is possible. Further calculations of the particle-hole self-energy are required to characterize this magnetic order. The competition between two types of susceptibility peaks mediating the d- and s-wave orders, shares many similarities with work done on iron-pnictide superconductors [Graser2009, Maier2011]. In iron pnictides a different nesting vector $\boldsymbol{Q}=(\pi,0)$ connects pockets and mediates the $s_{\pm}$ -wave. Ref. [Graser2009] predicts nearly degenerate multi-band $s$ - and $d$ -wave orders where a small change in the interaction parameters determines the favorable order. However, in contrast to our work their model does not contain spin-orbit coupling. In our work, the strong SOC and the consequent pseudospin degrees of freedom are partially responsible for the dominance of the $d$ -wave order in the realistic regime phase diagram.

The strain-induced superconducting regimes we find all occur when the Fermi surface has multiple pockets. All superconducting orders are thus multigap orders [Kresin2013, Zhitomirsky2001]. However, for a higher SOC or lower $U$ that would not necessarily be the case. The size of the second pocket FS₂ and the value of the Hund’s coupling determine the relative sizes of the gap for the two pockets. In Fig. 8, the relative gaps are shown projected onto the FS at low Hund’s coupling, and we find that pocket 1 dominates: Max $\left[\Delta_{\textrm{FS}_{1}}(\boldsymbol{k})\right]\gg$ Max $\left[\Delta_{\textrm{FS}_{2}}(\boldsymbol{k})\right]$ . The smaller gap is expected to have a smaller effect on the (shared) critical temperature. Methods to determine the relative size and pairing structures of two gap superconductors have been explored in multiple compounds such as MgB₂ [Souma2003] and SrTiO₃ [Fernandes2013]. Further proposals have been made to detect any offset phases between pairing functions. Especially for the $s_{\pm}$ -wave order where signatures of the phase difference between the two pockets could be detected via Josephson tunneling [Chen2010, Berg2011, Rodriguez-Mota2016].

Our results indicate that the doped and strained regime is of interest for potential iridate superconductivity. Understanding changes to magnetic fluctuations for any experiment combining strain and doping would provide great insight to the interplay of interactions and spin-orbit coupling in transition metal oxides. An observation of the susceptibility peaks that we identify as mediating superconductivity could hint at a possible superconducting phase nearby. Further understanding the signatures of the possible superconductivity, such as the Knight shift in the magnetic susceptibility [Romer2019] could be a direction for future work.

VII Acknowledgments

We acknowledge financial support from NSERC, RQMP, FRQNT, an Alexander McFee Fellowship from McGill University (LE), a Team Research Project from FRQNT, a grant from Fondation Courtois, and a Canada Research Chair. Computations were made on the supercomputers managed by Calcul Québec and Compute Canada. The operation of these supercomputers is funded by the Canada Foundation for Innovation (CFI), the ministère de l’économie, de la science et de l’innovation du Québec (MESI) and the Fonds de recherche du Québec - Nature et technologies (FRQ-NT).

Appendix A Phase diagrams at other U

In the model used in this paper there is a trade-off between doping- and strain- dependent behavior for a given strength of the Hubbard interaction $U$ , as motivated in section III. The value for the calculated the phase diagram, $U\approx 3|t|$ , finds a magnetic phase transition for an expected range of strain values. In Fig. 12 phase diagrams are calculated at a lower $U=0.78$ eV $\approx 2.2|t|$ . For the most realistic choice of Hund’s ( $J_{\textrm{H}}=0.1U$ ) and spin-orbit ( $\lambda=0.6$ eV), the magnetic region at $\epsilon=0$ only extends up to $n=4.85$ on the hole doped side. We note that for the chosen temperature and interaction parameters, a magnetic region extends beyond $n=5.3$ . Mean field studies for the lattice with staggered rotations, such as in Appendix B, reveal a possible canted ferromagnetic order for this doping. The high electron doping region is therefore likely to be a canted ferromagnet, as it has nesting vector $\boldsymbol{Q}=(\pi,\pi)$ . The phase diagram has a clear doping-dependence of the superconducting region, with a dome centered at the electron doped regime. However, the considered strain values only extend up to $\epsilon=-1.5\%$ .

For the choice $J_{\textrm{H}}=0.25U$ and $\lambda=0.5$ eV, neither a magnetic nor superconducting region is found for the temperature $T\approx 11$ K. Similarly to the result at $U\approx 3|t|$ , the leading eigenvector corresponds to two different types of pairing, depending on the considered doping regime.

To explore the robustness of the strain-induces $d$ -wave superconductivity, an additional phase diagram at $U\approx 4|t|$ is shown in Fig. 13. Even though the magnetic order persists up to higher values of the strain, the order is still a $d$ -wave with mainly $j=1/2$ contributions. At higher strains, the pocket FS₂ is larger than in Fig. 8 and the ratio Max $\left[\Delta_{\textrm{FS}_{2}}(\boldsymbol{k})\right]/$ Max $\left[\Delta_{\textrm{FS}_{1}}(\boldsymbol{k})\right]$ is increased.

Appendix B Mean field calculation

The phase diagrams presented in this work are limited in the determination of competition between orders. Even though all superconducting orders are treated equally, they are calculated in the normal state. Therefore any influence of particle-hole self-energies are ignored. To compare phase boundaries for the regions of interest a self-consistent mean-field calculation was performed. We have chosen to only compare the RPA calculation for the most realistic phase diagram in Fig. 3. For the $j=1/2$ superconducting order we can make an approximation of the effectively attractive interaction between sites. However, the anisotropic $s$ - and $d$ -wave orders, as well as for the $p$ -wave, require approximations of the effective attractive interactions between both nearest and next-nearest neighbor sites. As these multi-orbital pairing functions affect more than one $j$ -state, one must also approximate the strength of the interaction in each of these channels. The mean field calculation was chosen to include all possible on-site magnetic order parameters, in a two-site unit cell, as well as for a $d$ -wave superconducting order parameter in the $j=1/2$ state. The on-site order parameters, at sublattice $s=A,B$ , are

		$\displaystyle\langle c^{\dagger}_{s,a}c_{s,b}\rangle_{\textrm{MF},T}=$		(43)
		$\displaystyle\frac{1}{N_{k}^{2}}\sum_{n,\boldsymbol{k}}f(\xi_{\boldsymbol{k},n},T)\left[U_{\textrm{MF},\boldsymbol{k}}\right]_{(s,a)n}\left[U^{\dagger}_{\textrm{MF},\boldsymbol{k}}\right]_{n(s,b)}$		(43)

with the orbit-spin label $a=(\alpha,\sigma)$ , were calculated from a mean-field decoupling of the bare interactions in Eq. (20), as described in Ref. [25].

The superconducting order parameter is introduced by mean field decoupling of the Bogoliubov-de Genne (BdG) Hamiltonian. $\Delta_{j=1/2}$ is set to be a $j=1/2$ $d$ -wave singlet between the sublattices as

$\displaystyle H_{\text{MF,SC}}=\sum_{\boldsymbol{k}}$	$\displaystyle V(\epsilon)\Delta_{j=1/2}e^{ik_{x}}\left(\cos k_{x}-\cos k_{y}\right)$	(44)
	$\displaystyle\times\left[a_{\boldsymbol{k},A,(\frac{1}{2},+\frac{1}{2})}a_{-\boldsymbol{k},B,(\frac{1}{2},-\frac{1}{2})}\right.$
	$\displaystyle\left.-a_{\boldsymbol{k},A,(\frac{1}{2},-\frac{1}{2})}a_{-\boldsymbol{k},B,(\frac{1}{2},+\frac{1}{2})}+\text{h.c.}\right]$

where the operators $a$ are in the $j$ -state basis, as in section II, and $V(\epsilon)=-\frac{3}{4}J_{\text{eff}}(\epsilon)$ is the effective interaction [Laughlin2014, Farrell2014]. Due to the structure of the hopping terms we approximate $J_{\text{eff}}(\epsilon)$ in the $j=1/2$ subspace as $J_{\text{eff}}(\epsilon)=\sqrt{J_{1}^{2}+D^{2}}$ [31], where

J_{1}=\frac{4\left(t_{\text{eff}}(\epsilon)-t_{\text{eff},z}(\epsilon)\right)^{2}}{U_{\text{eff}}},\qquad D=\frac{8t_{\text{eff}}(\epsilon)t_{\text{eff},z}(\epsilon)}{U_{\text{eff}}}

(45)

Here the effective hopping parameters are $t_{\text{eff}}(\epsilon)=\frac{1}{3}\left(t_{1}+t_{4}+t_{5}\right)$ and $t_{\text{eff},z}(\epsilon)=t_{rot}$ . The order parameter is calculated via the self-consistency equation

	$\displaystyle\Delta_{j=1/2}=\frac{1}{N_{k}^{2}}\sum_{n,\boldsymbol{k}}\sum_{\alpha,\beta}$	$\displaystyle e^{-ik_{x}}M_{\left(A,(\frac{1}{2},+\frac{1}{2})\right),\alpha}M_{\left(B,(\frac{1}{2},-\frac{1}{2})\right),\beta}$		(46)
		$\displaystyle\times\left[U_{\textrm{MF},\boldsymbol{k}}\right]_{\alpha n}\left[U_{\textrm{MF},\boldsymbol{k}}\right]_{\beta n}f(\xi_{\boldsymbol{k},n},T)$		(46)

where the matrix $M$ is Eq. (8). This and Eq. (43) are solved simultaneously via iterations as a set of coupled equations. Here we are treating only the interaction within the $j=1/2$ as if it was a one-band model, where the interaction projected onto that subspace is $H_{\text{eff}}=\left(U-\frac{4}{3}J_{\text{H}}\right)n_{i,(\frac{1}{2},+\frac{1}{2})}n_{i,(\frac{1}{2},-\frac{1}{2})}=U_{\text{eff}}n_{i,(\frac{1}{2},+\frac{1}{2})}n_{i,(\frac{1}{2},-\frac{1}{2})}$ , so for the calculations $U_{\text{eff}}=\frac{13}{15}U$ . An inclusion of $d$ -wave pairing within the $(j,j_{z})=(\frac{3}{2},\pm\frac{3}{2})$ sate or between the $j=3/2$ and $j=1/2$ does not extend the superconducting phase in the calculation. The interaction used here is $U\approx t_{\text{eff}}$ , and not close to the strong coupling limit. However, to compare the competition between the $d$ -wave and magnetic orders for for the RPA calculations results in an order of roughly equal size.

Calculations were performed on a $N_{k}\times N_{k}$ lattice in momentum space, with $N_{k}=100$ , and self-consistent solutions were found iteratively with a convergence criterion of $10^{-7}$ . Due to the full set of order parameters containing terms that renormalize the spin-orbit coupling, the self-consistent mean field calculation predicts effects not included in the RPA calculation. However, the included mean field orders can only have nesting vectors $\boldsymbol{Q}=(0,0)$ or $\boldsymbol{Q}=(\pi,\pi)$ as the two sites allow us to study either net $(\langle c^{\dagger}_{A,a}c_{A,b}\rangle_{\textrm{MF},T}+\langle c^{\dagger}_{B,a}c_{B,b}\rangle_{\textrm{MF},T})/2$ or staggered $(\langle c^{\dagger}_{A,a}c_{A,b}\rangle_{\textrm{MF},T}-\langle c^{\dagger}_{B,a}c_{B,b}\rangle_{\textrm{MF},T})/2$ values of the order parameters.

In Fig. 14 there are two magnetic mean field orders. The most common order is the canted in-plane antiferromagnet (x-cAFM), with a staggered magnetic moment along the $x$ -direction and a net magnetic moment along the $y$ -direction. Similarly to previous works, the canting angle follows the rotations of the octahedra in the lattice. When a compressive strain increases the rotations, the magnetic canting angle thus increases as well. For a high enough $U$ the magnetic order remains up to high electron number ( $n>5.2$ ) the in-plane ferromagnetic component becomes dominant. A small AFM component remains, identifying this order as y-cFM. The canted orders consist mainly of $j=1/2$ pseudospins. However, as the SOC is lowered the AFM order in the hole doped region has small $(j,j_{z})=(\frac{3}{2},\pm\frac{3}{2})$ contributions. The order parameters which renormalize the SOC depend on strain and doping. However, if $\lambda=0.45$ eV is chosen the effective spin orbit coupling $\lambda_{\text{eff}}\approx 0.6$ eV.

A superconducting region is present in Fig. 14 and the competition between superconductivity and magnetism therefore does not affect its existence. However, in contrast to Fig. 3 superconductivity is only present for electron doping $n>5.05$ . For the non-interacting band structure, used for the RPA calculations in this paper, the $j=3/2$ pocket is present in the full charge doping region around $\epsilon=-2.2\%$ . As the $j=3/2$ hole pocket increases in size for a given charge doping, the $j=1/2$ electron pocket grows as well. The superconductivity in Fig. 3 is therefore present for a region where the $j=1/2$ electron pocket is significantly larger than for $\epsilon=0$ . There the $d$ -wave is present for all considered $n$ .

The discrepancy between the RPA and the mean field result is due to several factors. Since these two approaches make different approximations it is not possible to determine which phase diagram is more realistic. Instead, we can gain confidence in our results in parts of the phase diagram where the two approaches agree. Each approach has its strength and weaknesses. In mean field we must pre-determine the possible channels of superconductivity and the effective attractive interaction which does not change with strain and doping. On the other hand, the self-consistency equation Eq. (46) is not linearized like the RPA gap equation in Eq. (27) and therefore the mean field is better suited for determining the relative strength of the order parameters considered. The comparison between RPA and mean field theory therefore suggest that the electron doped region is more likely to host a d-wave superconducting order.

Appendix C Tetragonal splitting

Compression has an additional effect relevant to the iridates: an increased tetragonal distortion. The tetragonal distortion of the oxygen octahedra encompassing the iridium atoms has been measured to increase with compression - becoming more elongated when the in-plane compression increases. For unstrained Sr₂IrO₄ a small elongation is observed, which theoretically should result in a tetragonal splitting $\mu_{xy}>0$ . However, early ab initio calculations found that the band structure is best described by a shift $\mu_{xy}<0$ [6, Zhang2013]. Later works have proposed that the sign could arise from hybridization with ligand oxygen orbitals [Agrestini2017, Lane2020]. Due to this sign difference, previous works modeling the tetragonal splitting’s dependence on strain come to contradictory results, where $|\mu_{xy}|$ either increases or decreases [21, Zhang2013, 30]. A fitting of the change in energy splitting to RIXS measurements, in Ref. [28], found an increasing $|\mu_{xy}|$ for low compressive strain. A linearization of these results gives:

\mu_{xy}^{I}(\epsilon)=-|\mu_{xy}(0)|(1-0.2041\epsilon).

(47)

where $\epsilon$ is given in units of %. In Fig. 15 the largest eigenvalue for the linearized gap equation is compared for an approximation of $\mu_{xy}(\epsilon)$ as given by either Ref. [28] or Ref. [30]. The second approximation is based on theoretical calculations and the linearization is instead

\mu_{xy}^{D}(\epsilon)=-|\mu_{xy}(0)|(1+0.357\epsilon)

(48)

which results in a decreasing $|\mu_{xy}|$ under strain. We can note that even though the experimentally approximated values are only based on data points up to $\epsilon=-0.7$ %, the theoretical approximation $\mu_{xy}^{D}(\epsilon)$ is not compatible with the found trend. The experimentally motivated $\mu_{xy}^{I}(\epsilon)$ results in a slightly lower eigenvalue than the constant $\mu_{xy}$ for the $d$ -wave order, at $\lambda=0.6$ eV and $J_{\textrm{H}}=0.1U$ . Any change to the tetragonal splitting is thus expected to have a small impact on the strain-induced superconducting regions with $d$ -wave symmetry.

As a check of the magnetic phase boundaries, for different values of the tetragonal splitting, the mean field calculation in Appendix A was performed, for only the magnetic order parameters, at $N_{k}=200$ . In Fig. 16, an additional out-of-plane ferromagnetic order (z-FM) appears at the lower SOC and higher Hund’s coupling. This order is only favored at high enough $U$ and has a clear spin character, $\langle L_{z}\rangle/\langle S_{z}\rangle\approx 0.16$ . Contributions from each orbital are of approximately equal strength. In terms of $j$ -states the contributions are thus mainly from the $(j,j_{z})=(\frac{1}{2},\pm\frac{1}{2})$ and $(j,j_{z})=(\frac{3}{2},\pm\frac{3}{2})$ states.

As seen in Figs. 16 & 17, a larger absolute value of the splitting, $\mu_{xy}^{I}(\epsilon)$ , does not result in major changes to the magnetic phase boundaries. For the other option, $\mu_{xy}^{D}(\epsilon)$ , the behavior with doping changes as $\mu_{xy}^{D}(\epsilon\approx-2.8\%)=0$ . The x-cAFM order remains up to higher strains for hole instead of electron doping once the sign of $\mu_{xy}$ changes.

Appendix D Symmetry classification

In section V.1 the $\mathcal{S}\mathcal{P}\mathcal{O}\mathcal{T}$ -formalism [Sigrist1991, Linder2019] of classifying the symmetry of the superconducting pairing is introduced. The pairing must be antisymmetric under the product of operators

\mathcal{S}\mathcal{P}\mathcal{O}\mathcal{T}\Delta_{ab}(k)=\Delta_{ba}(-k)

(49)

where each exchange operator is defined as

$\displaystyle\mathcal{S}\Delta_{(s,\alpha,\sigma)(s^{\prime},\alpha^{\prime},\sigma^{\prime})}(\boldsymbol{k},\omega)=$	$\displaystyle\Delta_{(s,\alpha,\sigma^{\prime})(s^{\prime},\alpha^{\prime},\sigma)}(\boldsymbol{k},\omega)$	(50)
$\displaystyle\mathcal{P}\Delta_{(s,\alpha,\sigma)(s^{\prime},\alpha^{\prime},\sigma^{\prime})}(\boldsymbol{k},\omega)=$	$\displaystyle\Delta_{(s,\alpha,\sigma)(s^{\prime},\alpha^{\prime},\sigma^{\prime})}(-\boldsymbol{k},\omega)$
$\displaystyle\mathcal{O}\Delta_{(s,\alpha,\sigma)(s^{\prime},\alpha^{\prime},\sigma^{\prime})}(\boldsymbol{k},\omega)=$	$\displaystyle\Delta_{(s,\alpha^{\prime},\sigma)(s^{\prime},\alpha,\sigma^{\prime})}(\boldsymbol{k},\omega)$
$\displaystyle\mathcal{T}\Delta_{(s,\alpha,\sigma)(s^{\prime},\alpha^{\prime},\sigma^{\prime})}(\boldsymbol{k},\omega)=$	$\displaystyle\Delta_{(s,\alpha,\sigma)(s^{\prime},\alpha^{\prime},\sigma^{\prime})}(\boldsymbol{k},-\omega)$

In Fig. 18 the symmetries of a found $j=1/2$ $d$ -wave pairing is shown. The maxima of each component of $\Delta$ are separated under the present spin-singlet and spin-triplet operations. An ideal pseudospin singlet will have both spin-singlet and spin-triplet components of comparable size, in the orbital-spin basis. As the transformation between bases is independent of momentum, the spatial parity of the state is unchanged. The spin-triplet components are therefore orbital-singlets $\mathcal{S}^{+}\mathcal{P}^{+}\mathcal{O}^{-}\mathcal{T}^{+}$ . The found $d$ -wave in Fig. 18 has additional small components from other $j$ -states. The compressive strain decreases the contributions from the $xy$ -orbital at the Fermi surface as well as increases the interorbital $xz$ - $yz$ hopping. There are therefore more contributions from spin-singlet and $m_{z}=0$ spin-triplet components, than from $m_{z}=1$ triplets. In Fig. 19 the symmetries of an $s_{\pm}$ -wave is shown. The $s_{\pm}$ -wave has strong components from several $j$ -states. As this order is found to be mediated by spin-like fluctuations mainly in the $yz$ - and $xz$ -orbitals, between bands of $j=1/2$ and $j=3/2$ character, we can consider the main components of the pairing to come from pairing within and between those orbitals.

The superconducting order can be expressed exactly via its full symmetry representation. We consider the pairing for a spin and orbital combination:

\Delta_{\sigma\sigma^{\prime}}^{\alpha\beta}(\boldsymbol{k})=\frac{1}{2}\sum_{s,s^{\prime}}e^{ik_{x}(\Theta_{s}-\Theta_{s^{\prime}})}\Delta_{\sigma\sigma^{\prime}}^{(\alpha,s)(\beta,s^{\prime})}(\boldsymbol{k})

(51)

where $\Theta_{s}$ is the same function as in Eq. (38) and gives us the combined contribution from both sublattice sites. The symmetry representation for the spatial symmetry $\eta^{\mu}(\boldsymbol{k})$ is specified in Table 1. The pairing can be decomposed into symmetry representations for the spin and orbital structure. If $\mathcal{C}_{\mu\nu\rho}$ is the projection constant for a chosen representation, then any superconducting order can be written as:

\Delta_{\sigma\sigma^{\prime}}^{\alpha\beta}(\boldsymbol{k})=\sum_{\mu,\nu,\rho}\mathcal{C}_{\mu\nu\rho}\eta^{\mu}(\boldsymbol{k})S^{\nu}_{\sigma\sigma^{\prime}}O^{\rho}_{\alpha\beta}

(52)

For the spin degree of freedom, $S^{\nu}_{\sigma\sigma^{\prime}}$ are the generators for the SU(2) algebra, the Pauli matrices $\sigma^{i}$ with $i=0,x,y,z$ . For the orbital degree of freedom, $O^{\rho}_{\alpha\beta}$ are the generators for the SU(3) algebra, the Gell-Mann matrices [GellMann1961] $\lambda_{i}$ with $i=0,1,\dots,8$ . The Pauli matrices act in spin $(\uparrow,\downarrow)$ space and can form spin-triplets $(\sigma^{0},\sigma^{x},\sigma^{z})$ and spin-singlets $(\sigma^{y})$ . The Gell-Mann matrices act in orbital $(d_{yz},d_{xz},d_{xy})$ space and can form intraorbital pairings $(\lambda_{0},\lambda_{3},\lambda_{8})$ , and interorbital pairings that can be either even $(\lambda_{1},\lambda_{4},\lambda_{6})$ or odd $(\lambda_{2},\lambda_{5},\lambda_{7})$ under orbital exchange.

The found $d$ -wave is a pseudospin singlet, which within the $j=1/2$ subspace is

\Delta_{d}(\boldsymbol{k})\approx\eta^{B_{1g}}_{R=1}(\boldsymbol{k})\otimes(i\tilde{\sigma}^{y})

(53)

where $\tilde{\sigma}^{y}$ acts on the pseudospins $\tau=+,-$ . The $s_{\pm}$ -wave can be expressed approximately as intraorbital spin-singlets and interorbital spin-triplets in the $yz$ - and $xz$ -orbitals:

	$\displaystyle\Delta_{s_{\pm}}(\boldsymbol{k})\approx$	$\displaystyle\left(\mathcal{C}^{A_{1g}}_{R=0}\eta^{A_{1g}}_{R=0}(\boldsymbol{k})+\eta^{A_{1g}}_{R=2}(\boldsymbol{k})\right)\otimes$		(54)
		$\displaystyle\left[(i\sigma^{y})\otimes\left(\lambda_{11}+\lambda_{22}\right)-i\sigma^{x}\otimes\left(i\lambda_{2}\right)\right]$		(54)

with $\lambda_{11}=\left(\frac{1}{3}\lambda_{0}+\frac{1}{2}\lambda_{3}+\frac{1}{2\sqrt{3}}\lambda_{8}\right)$ and $\lambda_{22}=\left(\frac{1}{3}\lambda_{0}-\frac{1}{2}\lambda_{3}+\frac{1}{2\sqrt{3}}\lambda_{8}\right)$ representing intraorbital pairing. The A_1g spatial symmetry has two contributions with a relative weight specified by $\mathcal{C}^{A_{1g}}_{R=0}\approx 0.7$ . This pairing has additional smaller contributions involving the $xy$ -orbital: $\propto(i\sigma^{y})\otimes\lambda_{33}=(i\sigma^{y})\otimes\left(\lambda_{0}-\sqrt{3}\lambda_{8}\right)$ , $\propto(-i)\sigma^{0}\otimes\left(i\lambda_{5}\right)$ , and $\propto(-i)\sigma^{z}\otimes\left(i\lambda_{7}\right)$ .

Appendix E Anisotropic pairing

The anisotropic orders found at high compressive strain, in section V.3, appear when fluctuations with $\boldsymbol{Q}$ -vectors which connect states of the same orbitals become large enough. In Fig. 20 the found anisotropic orders, for three values of the Hund’s coupling, are characterized by projecting the pairing on the Fermi surface as well as the maximal value of each orbital-spin component of the pairing. At the Fermi surface the weight is stronger along the $x$ -direction, on all pockets. The $xz$ orbital has the largest contribution. The other notable difference, from the $d$ - and $s_{\pm}$ -wave orders, is that the anisotropic orders have equal size intra- and inter-sublattice components of the $xz$ -orbital.

The anisotropic orders and $s_{\pm}$ -order are mediated mainly by spin fluctuations. However, the pockets at the Fermi surface have a strong character of the $(j,j_{z})=(\frac{1}{2},\pm\frac{1}{2})$ or $(j,j_{z})=(\frac{3}{2},\pm\frac{3}{2})$ state. The orders originate mostly from the $yz$ - and $xz$ -orbitals, and the pairing can be transformed into the $j$ -states, via Eq. (8), as

	$\displaystyle\Delta^{\tau\tau^{\prime}}_{mn}(k)=\sum_{\alpha\beta\sigma\sigma^{\prime}}c^{(\tau\tau^{\prime})(\sigma\sigma^{\prime})}_{(mn)(\alpha\beta)}\Delta^{\sigma\sigma^{\prime}}_{\alpha\beta}(k)$		(55)
	$\displaystyle=\sum_{\alpha\beta\sigma\sigma^{\prime}}M_{(m,\tau),(\alpha,\sigma)}M_{(n,\tau^{\prime}),(\beta,\sigma^{\prime})}\Delta^{\sigma\sigma^{\prime}}_{\alpha\beta}(k)$		(55)

where $m=\left(\frac{1}{2},\pm\frac{1}{2}\right),\left(\frac{3}{2},\pm\frac{1}{2}\right),\left(\frac{3}{2},\pm\frac{3}{2}\right)$ , $\tau=+,-$ , $\alpha=yz,xz,xy$ , and $\sigma=\uparrow,\downarrow$ . We can identify $c^{(+-)(\uparrow\downarrow)}_{(\frac{3}{2},\pm\frac{3}{2})(yz,yz)}=c^{(+-)(\uparrow\downarrow)}_{(\frac{3}{2},\pm\frac{3}{2})(xz,xz)}=-\frac{1}{2}$ and $c^{(+-)(\uparrow\downarrow)}_{(\frac{1}{2},\pm\frac{1}{2})(yz,yz)}=c^{(+-)(\uparrow\downarrow)}_{(\frac{1}{2},\pm\frac{1}{2})(xz,xz)}=\frac{1}{3}$ . As can be seen in Fig. 5, sections of the hole doped Fermi surface are mainly of either $yz$ - or $xz$ -character. This is a result of the quasi-1d dispersion in each of these orbitals in Eq. (4). A simplified model for the two Fermi surfaces is to introduce orbitals with a spatial dependence along the parameter $\theta$ around a circular FS:

$\displaystyle\|yz\rangle_{\theta}=$	$\displaystyle\|\cos\theta\|\|yz\rangle$	(56)
$\displaystyle\|xz\rangle_{\theta}=$	$\displaystyle\|\sin\theta\|\|xz\rangle$
$\displaystyle\|xy\rangle_{\theta}=$	$\displaystyle\left(\|\cos\theta\|+\|\sin\theta\|\right)\|xy\rangle$

Here only one site and the full Brillouin zone (BZ) are considered for simplicity. So in a BZ where FS₁ is purely $(\frac{1}{2},\pm\frac{1}{2})$ and FS₂ is $(\frac{3}{2},\pm\frac{3}{2})$ (with an energy shift $\xi$ of the $xy$ -orbital):

	$\displaystyle\|\textrm{FS}_{1},\uparrow\rangle_{\theta}=\|\frac{1}{2},+\frac{1}{2}\rangle_{\theta}-\xi\|xy,\uparrow\rangle_{\theta}$		(57)
	$\displaystyle=\frac{1}{\mathcal{N}_{1}}\left(\|yz,\downarrow\rangle_{\theta}-i\|xz,\downarrow\rangle_{\theta}+(1-\sqrt{3}\xi)\|xy,\uparrow\rangle_{\theta}\right)$		(57)

and

|\textrm{FS}_{2},\uparrow\rangle_{\theta}=|\frac{3}{2},+\frac{3}{2}\rangle_{\theta}=\frac{1}{\mathcal{N}_{2}}\left(-|yz,\downarrow\rangle_{\theta}-i|xz,\downarrow\rangle_{\theta}\right)

(58)

where $\mathcal{N}_{1},\mathcal{N}_{2}$ are normalization factors. For a simplified constant order only within each orbital, $\Delta^{\uparrow\downarrow}_{yz,yz}(k)=\Delta$ the pairing on each FS becomes

\Delta^{\downarrow\uparrow}_{\textrm{FS}_{1}}=\frac{\cos^{2}\theta}{\mathcal{N}_{1}^{2}}\Delta,\qquad\Delta^{\downarrow\uparrow}_{\textrm{FS}_{2}}=-\frac{\cos^{2}\theta}{\mathcal{N}_{2}^{2}}\Delta

(59)

If instead $\Delta^{\uparrow\downarrow}_{xz,xz}=\Delta$

\Delta^{\downarrow\uparrow}_{\textrm{FS}_{1}}=\frac{\sin^{2}\theta}{\mathcal{N}_{1}^{2}}\Delta,\qquad\Delta^{\downarrow\uparrow}_{\textrm{FS}_{2}}=-\frac{\sin^{2}\theta}{\mathcal{N}_{2}^{2}}\Delta

(60)

Each orbital thus results in a quasi-1d pairing on both Fermi surfaces. If one of the orbitals dominate $\Delta^{\uparrow\downarrow}_{xz,xz}>\Delta^{\uparrow\downarrow}_{yz,yz}$ the order is an anisotropic $s$ -wave, with some $d$ -wave components.

To study the full $s_{\pm}$ -wave, it can be modeled as equal parts from both orbitals, $\Delta^{\uparrow\downarrow}_{yz,yz}=\Delta^{\uparrow\downarrow}_{xz,xz}=\Delta$ :

	$\displaystyle\Delta^{\downarrow\uparrow}_{\textrm{FS}_{1}}=$	$\displaystyle\frac{\cos^{2}\theta+\sin^{2}\theta}{\mathcal{N}_{1}^{2}}\Delta$	$\displaystyle=\frac{1}{\mathcal{N}_{1}^{2}}\Delta$		(61)
	$\displaystyle\Delta^{\downarrow\uparrow}_{\textrm{FS}_{2}}=$	$\displaystyle-\frac{\cos^{2}\theta+\sin^{2}\theta}{\mathcal{N}_{2}^{2}}\Delta$	$\displaystyle=-\frac{1}{\mathcal{N}_{2}^{2}}\Delta$		(61)

where in the ideal $j$ -state case $\mathcal{N}_{1}^{2}=3$ and $\mathcal{N}_{2}^{2}=2$ . This is one of the primary reasons for the $s$ -wave symmetry resulting in an $s_{\pm}$ -wave with a larger weight on FS₂. However, the found $s_{\pm}$ -wave does not have purely $\mathcal{S}^{-}\mathcal{P}^{+}\mathcal{O}^{+}\mathcal{T}^{+}$ intraorbital components but also interorbital $\mathcal{S}^{+}\mathcal{P}^{+}\mathcal{O}^{-}\mathcal{T}^{+}$ -terms, see Fig. 19. As these terms all have the same magnitude $\Delta^{\uparrow\downarrow}_{yz,xz}=-\Delta^{\uparrow\downarrow}_{xz,yz}=i\Delta$ . Projected onto the Fermi surface

\Delta^{\downarrow\uparrow}_{\textrm{FS}_{1}}=\frac{2|\cos\theta||\sin\theta|}{\mathcal{N}_{1}^{2}}\Delta,\qquad\Delta^{\downarrow\uparrow}_{\textrm{FS}_{2}}=\frac{2|\cos\theta||\sin\theta|}{\mathcal{N}_{2}^{2}}\Delta.

(62)

They contribute to both bands with the same sign and to the same sections. The pairing used for these examples so far has been a uniform $s$ -wave, $\Delta(\boldsymbol{k})=\Delta$ . In Eq. (54) we can note that the found $s_{\pm}$ -wave has a dependence on momentum, with a large contribution from the $\eta^{A_{1g}}_{R=2}(\boldsymbol{k})=2\cos k_{x}\cos k_{y}$ symmetry. The placement of FS₁ and FS₂ in the BZ therefore affects the sign of the gaps. As a result, both intra- and inter-orbital terms play a role in the origin of the $s_{\pm}$ -wave pairing.

Appendix F Topological invariant

To determine the topological properties of the found odd parity order we consider the pairing in the Bogoliubov-de Gennes (BdG) Hamiltonian

H_{\textrm{BdG}}=\sum_{\boldsymbol{k}}\left(c^{\dagger}_{\boldsymbol{k}},c_{-\boldsymbol{k}}\right)\left(\begin{array}[]{cc}H(\boldsymbol{k})&\Delta(\boldsymbol{k})\\ \Delta^{\dagger}(\boldsymbol{k})&-H^{\textrm{T}}(-\boldsymbol{k})\end{array}\right)\left(\begin{array}[]{c}c_{\boldsymbol{k}}\\ c^{\dagger}_{-\boldsymbol{k}}\end{array}\right)

(63)

with the non-interacting Hamiltonian $H(\boldsymbol{k})$ from Eq. (1). $\Delta(\boldsymbol{k})$ is the eigenstate of the largest eigenvalue for Eq. (29), scaled to set the minimal gap at the Fermi surface to 0.02eV. As the specific odd parity pairing found in section V.5 is block diagonal in pseudospin up and down, the full Hamiltonian can be rearranged into a block-diagonal form

H_{\textrm{BdG}}=\sum_{\boldsymbol{k}}\boldsymbol{\Psi}_{\boldsymbol{k}}^{\dagger}\left(\begin{array}[]{cc}H_{\textrm{BdG}}^{+}(\boldsymbol{k})&0\\ 0&H_{\textrm{BdG}}^{-}(\boldsymbol{k})\end{array}\right)\boldsymbol{\Psi}_{\boldsymbol{k}}

(64)

where $\boldsymbol{\Psi}_{\boldsymbol{k}}=\left(\Psi_{+,\boldsymbol{k}},\Psi_{-,\boldsymbol{k}}\right)^{\textrm{T}}$ is divided into the pseudospin sectors $\{(yz,\downarrow),(xz,\downarrow),(xy,\uparrow)\}$ and $\{(yz,\uparrow),(xz,\uparrow),(xy,\downarrow)\}$ . The eigenstates and eigenvalues are given as $H_{\textrm{BdG}}|n\rangle=E_{n}|n\rangle$ and for each pseudospin sector $H_{\textrm{BdG}}^{\tau}|n\rangle=E_{n,\tau}|n,\tau\rangle$ . We calculate the Berry curvature for all filled bands via [Berry1984]

	$\displaystyle\Omega^{(n)}_{z}(\boldsymbol{k})=$	$\displaystyle-\textrm{Im}\sum_{n^{\prime}\neq n}f(E_{n})f(1-E_{n^{\prime}})$		(65)
		$\displaystyle\times\frac{\langle n\|\partial_{x}H_{\textrm{BdG}}\|n^{\prime}\rangle\langle n^{\prime}\|\partial_{y}H_{\textrm{BdG}}\|n\rangle-\langle n\|\partial_{y}H_{\textrm{BdG}}\|n^{\prime}\rangle\langle n^{\prime}\|\partial_{x}H_{\textrm{BdG}}\|n\rangle}{(E_{n}-E_{n^{\prime}})^{2}}$		(65)

	$\displaystyle\Omega^{(n)}_{z}(\boldsymbol{k})=$	$\displaystyle-\textrm{Im}\sum_{n^{\prime}\neq n}\sum_{a,b,c,d}\sum_{i,j\in\{x,y\}}f(E_{n})f(1-E_{n^{\prime}})\left[U_{\boldsymbol{k}}^{\dagger}\right]_{na}\left[U_{\boldsymbol{k}}\right]_{bn^{\prime}}\left[U_{\boldsymbol{k}}^{\dagger}\right]_{n^{\prime}c}\left[U_{\boldsymbol{k}}\right]_{dn}$		(66)
		$\displaystyle\times\epsilon_{ij}\frac{\langle a\|\partial_{i}H_{\textrm{BdG}}\|b\rangle\langle c\|\partial_{j}H_{\textrm{BdG}}\|d\rangle}{(E_{n}-E_{n^{\prime}})^{2}}$		(66)

where $\epsilon_{xy}=-\epsilon_{yx}=1$ and $\epsilon_{ii}=0$ . $|a\rangle$ is the spin-orbital basis $a=(\alpha,\sigma,s)$ with orbital $\alpha$ , spin $\sigma$ , and sublattice $s$ . However, since all $\langle a,+|\partial_{i}H_{\textrm{BdG}}|b,-\rangle=0$ the Berry curvature can be calculated separately for each pseudospin sector

\Omega^{(n)}_{z}(\boldsymbol{k})=\Omega^{(n,+)}_{z}(\boldsymbol{k})+\Omega^{(n,-)}_{z}(\boldsymbol{k})

(67)

and the Chern number per pseudospin is

\mathcal{C}_{\tau}=\frac{1}{2\pi}\int_{\textrm{BZ}}\sum_{n}\Omega^{(n,\tau)}(\boldsymbol{k})dk_{x}dk_{y}

(68)

In Fig. 21 the Berry curvature for all filled bands is shown for the pseudospin down sector. In addition, the phase $\phi_{\boldsymbol{k}}$ of the intraband pairing, where $\Delta_{\textrm{FS}_{n}}^{\downarrow\downarrow}(\boldsymbol{k})=|\Delta_{\textrm{FS}_{n}}^{\downarrow\downarrow}(\boldsymbol{k})|e^{i\phi_{\boldsymbol{k}}}$ , is plotted along the Fermi surface. The winding of the phase can be seen to be opposite for the bands. If the pockets are fully separated the Chern number for each pocket can also be given by the winding of the phase [Cheng2010]

\mathcal{C}_{n}=\frac{1}{2\pi}\oint_{\textrm{FS}_{n}}\nabla\phi_{\boldsymbol{k}}\cdot d\boldsymbol{k}.

(69)

References

Witczak-Krempa et al. [2014] W. Witczak-Krempa, G. Chen, Y. B. Kim, and L. Balents, Correlated Quantum Phenomena in the Strong Spin-Orbit Regime, Annual Review of Condensed Matter Physics 5, 57 (2014).
Cao and Schlottmann [2018] G. Cao and P. Schlottmann, The challenge of spin–orbit-tuned ground states in iridates: a key issues review, Reports on Progress in Physics 81, 042502 (2018).
Bertinshaw et al. [2019] J. Bertinshaw, Y. Kim, G. Khaliullin, and B. Kim, Square Lattice Iridates, Annual Review of Condensed Matter Physics 10, 315 (2019).
Boseggia et al. [2013] S. Boseggia, H. C. Walker, J. Vale, R. Springell, Z. Feng, R. S. Perry, M. Moretti Sala, H. M. Rønnow, S. P. Collins, and D. F. McMorrow, Locking of iridium magnetic moments to the correlated rotation of oxygen octahedra in Sr₂IrO₄ revealed by x-ray resonant scattering, Journal of Physics: Condensed Matter 25, 422202 (2013).
Crawford et al. [1994] M. K. Crawford, M. A. Subramanian, R. L. Harlow, J. A. Fernandez-Baca, Z. R. Wang, and D. C. Johnston, Structural and magnetic studies of Sr₂IrO₄, Phys. Rev. B 49, 9198 (1994).
Kim et al. [2008] B. J. Kim, H. Jin, S. J. Moon, J.-Y. Kim, B.-G. Park, C. S. Leem, J. Yu, T. W. Noh, C. Kim, S.-J. Oh, J.-H. Park, V. Durairaj, G. Cao, and E. Rotenberg, Novel $J_{\text{eff}}=1/2$ Mott State Induced by Relativistic Spin-Orbit Coupling in Sr₂IrO₄, Physical Review Letters 101, 076402 (2008).
Kim et al. [2014] Y. K. Kim, O. Krupin, J. D. Denlinger, A. Bostwick, E. Rotenberg, Q. Zhao, J. F. Mitchell, J. W. Allen, and B. J. Kim, Fermi arcs in a doped pseudospin-1/2 Heisenberg antiferromagnet, Science 345, 187 (2014).
Cao et al. [2016] Y. Cao, Q. Wang, J. A. Waugh, T. J. Reber, H. Li, X. Zhou, S. Parham, S. R. Park, N. C. Plumb, E. Rotenberg, A. Bostwick, J. D. Denlinger, T. Qi, M. A. Hermele, G. Cao, and D. S. Dessau, Hallmarks of the Mott-metal crossover in the hole-doped pseudospin-1/2 Mott insulator Sr₂IrO₄, Nature Communications 7, 1 (2016).
Wang and Senthil [2011] F. Wang and T. Senthil, Twisted Hubbard Model for Sr₂IrO₄ : Magnetism and Possible High Temperature Superconductivity, Physical Review Letters 106, 136402 (2011).
Watanabe et al. [2013] H. Watanabe, T. Shirakawa, and S. Yunoki, Monte Carlo Study of an Unconventional Superconducting Phase in Iridium Oxide $J_{\text{eff}}=1/2$ Mott Insulators Induced by Carrier Doping, Physical Review Letters 110, 027002 (2013).
Yang et al. [2014] Y. Yang, W.-S. Wang, J.-g. Liu, H. Chen, J.-h. Dai, and Q.-h. Wang, Superconductivity in doped Sr₂IrO₄: A functional renormalization group study, Physical Review B 89, 094518 (2014).
Meng et al. [2014] Z. Y. Meng, Y. B. Kim, and H.-Y. Kee, Odd-Parity Triplet Superconducting Phase in Multiorbital Materials with a Strong Spin-Orbit Coupling: Application to Doped ${\mathrm{Sr}}_{2}{\mathrm{IrO}}_{4}$ , Phys. Rev. Lett. 113, 177003 (2014).
Gao et al. [2015] Y. Gao, T. Zhou, H. Huang, and Q. H. Wang, Possible superconductivity in Sr₂IrO₄ probed by quasiparticle interference, Scientific Reports 5, 1 (2015).
Nishiguchi et al. [2019] K. Nishiguchi, T. Shirakawa, H. Watanabe, R. Arita, and S. Yunoki, Possible Superconductivity Induced by a Large Spin–Orbit Coupling in Carrier Doped Iridium Oxide Insulators: A Weak Coupling Approach, Journal of the Physical Society of Japan 88, 094701 (2019).
Kim et al. [2012] J. Kim, D. Casa, M. H. Upton, T. Gog, Y.-J. Kim, J. F. Mitchell, M. van Veenendaal, M. Daghofer, J. van den Brink, G. Khaliullin, and B. J. Kim, Magnetic Excitation Spectra of Sr₂IrO₄ Probed by Resonant Inelastic X-Ray Scattering: Establishing Links to Cuprate Superconductors, Phys. Rev. Lett. 108, 177003 (2012).
Kim et al. [2016] Y. K. Kim, N. H. Sung, J. D. Denlinger, and B. J. Kim, Observation of a d-wave gap in electron-doped Sr₂IrO₄, Nature Physics 12, 37 (2016).
Chen et al. [2015] X. Chen, T. Hogan, D. Walkup, W. Zhou, M. Pokharel, M. Yao, W. Tian, T. Z. Ward, Y. Zhao, D. Parshall, C. Opeil, J. W. Lynn, V. Madhavan, and S. D. Wilson, Influence of electron doping on the ground state of $({\mathrm{Sr}}_{1-x}{\mathrm{La}}_{x}{)}_{2}{\mathrm{IrO}}_{4}$ , Phys. Rev. B 92, 075125 (2015).
Battisti et al. [2017] I. Battisti, K. M. Bastiaans, V. Fedoseev, A. de la Torre, N. Iliopoulos, A. Tamai, E. C. Hunter, R. S. Perry, J. Zaanen, F. Baumberger, and M. P. Allan, Universality of pseudogap and emergent order in lightly doped Mott insulators, Nature Physics 13, 21 (2017).
Chen et al. [2018] X. Chen, J. L. Schmehr, Z. Islam, Z. Porter, E. Zoghlin, K. Finkelstein, J. P. Ruff, and S. D. Wilson, Unidirectional spin density wave state in metallic (Sr_1-xLa_x)₂IrO₄, Nature Communications 9, 103 (2018).
Lindquist and Kee [2019] A. W. Lindquist and H.-Y. Kee, Odd-parity superconductivity driven by octahedra rotations in iridium oxides, Physical Review B 100, 054512 (2019).
Jackeli and Khaliullin [2009] G. Jackeli and G. Khaliullin, Mott Insulators in the Strong Spin-Orbit Coupling Limit: From Heisenberg to a Quantum Compass and Kitaev Models, Phys. Rev. Lett. 102, 017205 (2009).
Perkins et al. [2014] N. B. Perkins, Y. Sizyuk, and P. Wölfle, Interplay of many-body and single-particle interactions in iridates and rhodates, Physical Review B 89, 035143 (2014).
Liu et al. [2015] P. Liu, S. Khmelevskyi, B. Kim, M. Marsman, D. Li, X.-Q. Chen, D. D. Sarma, G. Kresse, and C. Franchini, Anisotropic magnetic couplings and structure-driven canted to collinear transitions in ${\mathrm{Sr}}_{2}{\mathrm{IrO}}_{4}$ by magnetically constrained noncollinear DFT, Phys. Rev. B 92, 054428 (2015).
Torchinsky et al. [2015] D. H. Torchinsky, H. Chu, L. Zhao, N. B. Perkins, Y. Sizyuk, T. Qi, G. Cao, and D. Hsieh, Structural Distortion-Induced Magnetoelastic Locking in Sr₂IrO₄ Revealed through Nonlinear Optical Harmonic Generation, Physical Review Letters 114, 096404 (2015).
Engström et al. [2021] L. Engström, T. Pereg-Barnea, and W. Witczak-Krempa, Modeling multiorbital effects in Sr₂IrO₄ under strain and a Zeeman field, Physical Review B 103, 155147 (2021).
Lupascu et al. [2014] A. Lupascu, J. P. Clancy, H. Gretarsson, Z. Nie, J. Nichols, J. Terzic, G. Cao, S. S. A. Seo, Z. Islam, M. H. Upton, J. Kim, D. Casa, T. Gog, A. H. Said, V. M. Katukuri, H. Stoll, L. Hozoi, J. van den Brink, and Y.-J. Kim, Tuning Magnetic Coupling in Sr₂IrO₄ Thin Films with Epitaxial Strain, Physical Review Letters 112, 147201 (2014).
Hao et al. [2019] L. Hao, D. Meyers, M. Dean, and J. Liu, Novel spin-orbit coupling driven emergent states in iridate-based heterostructures, Journal of Physics and Chemistry of Solids 128, 39 (2019).
Paris et al. [2020] E. Paris, Y. Tseng, E. M. Pärschke, W. Zhang, M. H. Upton, A. Efimenko, K. Rolfs, D. E. McNally, L. Maurel, M. Naamneh, M. Caputo, V. N. Strocov, Z. Wang, D. Casa, C. W. Schneider, E. Pomjakushina, K. Wohlfeld, M. Radovic, and T. Schmitt, Strain engineering of the charge and spin-orbital interactions in Sr₂IrO₄, Proceedings of the National Academy of Sciences 117, 24764 (2020).
Samanta et al. [2018] K. Samanta, F. M. Ardito, N. M. Souza-Neto, and E. Granado, First-order structural transition and pressure-induced lattice/phonon anomalies in Sr₂IrO₄, Physical Review B 98, 094101 (2018).
Bhandari et al. [2019] C. Bhandari, Z. S. Popović, and S. Satpathy, Electronic structure and optical properties of Sr₂IrO₄ under epitaxial strain, New Journal of Physics 21, 10.1088/1367-2630/aaff1b (2019).
Seo et al. [2019] A. Seo, P. P. Stavropoulos, H.-H. Kim, K. Fürsich, M. Souri, J. G. Connell, H. Gretarsson, M. Minola, H. Y. Kee, and B. Keimer, Compressive strain induced enhancement of exchange interaction and short-range magnetic order in $\mathrm{S}{\mathrm{r}}_{2}\mathrm{Ir}{\mathrm{O}}_{4}$ investigated by Raman spectroscopy, Phys. Rev. B 100, 165106 (2019).
Trommler et al. [2010] S. Trommler, R. Hühne, K. Iida, P. Pahlke, S. Haindl, L. Schultz, and B. Holzapfel, Reversible shift in the superconducting transition for La_1.85Sr_0.15CuO₄ and BaFe_1.8Co_0.2As₂ using piezoelectric substrates, New Journal of Physics 12, 0 (2010).
Hühne et al. [2008] R. Hühne, D. Okai, K. Dörr, S. Trommler, A. Herklotz, B. Holzapfel, and L. Schultz, Dynamic investigations on the influence of epitaxial strain on the superconducting transition in YBa₂Cu₃O_7-x, Superconductor Science and Technology 21, 075020 (2008).
Phan et al. [2017] G. N. Phan, K. Nakayama, K. Sugawara, T. Sato, T. Urata, Y. Tanabe, K. Tanigaki, F. Nabeshima, Y. Imai, A. Maeda, and T. Takahashi, Effects of strain on the electronic structure, superconductivity, and nematicity in FeSe studied by angle-resolved photoemission spectroscopy, Physical Review B 95, 224507 (2017).
Ahadi et al. [2019] K. Ahadi, L. Galletti, Y. Li, S. Salmani-Rezaie, W. Wu, and S. Stemmer, Enhancing superconductivity in SrTiO₃ films with strain, Science Advances 5, 1 (2019).
Zhang et al. [2022] J. Zhang, H. Wu, G. Zhao, L. Han, and J. Zhang, A review on strain study of cuprate superconductors, Nanomaterials 12, 10.3390/nano12193340 (2022).
Kawasugi et al. [2008] Y. Kawasugi, H. M. Yamamoto, M. Hosoda, N. Tajima, T. Fukunaga, K. Tsukagoshi, and R. Kato, Strain-induced superconductor/insulator transition and field effect in a thin single crystal of molecular conductor, Applied Physics Letters 92, 24 (2008).
Liu et al. [2018] Y.-C. Liu, W.-S. Wang, F.-C. Zhang, and Q.-H. Wang, Superconductivity in ${\mathrm{Sr}}_{2}{\mathrm{RuO}}_{4}$ thin films under biaxial strain, Phys. Rev. B 97, 224522 (2018).
Haskel et al. [2020] D. Haskel, G. Fabbris, J. H. Kim, L. S. I. Veiga, J. R. L. Mardegan, C. A. Escanhoela, S. Chikara, V. Struzhkin, T. Senthil, B. J. Kim, G. Cao, and J.-W. Kim, Possible Quantum Paramagnetism in Compressed Sr₂IrO₄, Physical Review Letters 124, 067201 (2020).
Chen et al. [2020] C. Chen, Y. Zhou, X. Chen, T. Han, C. An, Y. Zhou, Y. Yuan, B. Zhang, S. Wang, R. Zhang, L. Zhang, C. Zhang, Z. Yang, L. E. DeLong, and G. Cao, Persistent insulating state at megabar pressures in strongly spin-orbit coupled Sr₂IrO₄, Physical Review B 101, 144102 (2020).
Li et al. [2021]

$\displaystyle\|yz\rangle_{\theta}=$	$\displaystyle\|\cos\theta\|\|yz\rangle$	(56)
$\displaystyle\|xz\rangle_{\theta}=$	$\displaystyle\|\sin\theta\|\|xz\rangle$
$\displaystyle\|xy\rangle_{\theta}=$	$\displaystyle\left(\|\cos\theta\|+\|\sin\theta\|\right)\|xy\rangle$

Strain-induced superconductivity in Sr2IrO4