# Orbitally induced string formation in the spin-orbital polarons

##### I. INTRODUCTION

A single hole doped into the half-filled Mott insulating ground state with antiferromagnetic (AF) order does not move freely as its motion disturbs the spin order. ^{1} Instead, it couples to the collective (delocalized) excitations of the AF ordered phase (magnons), and it propagates through the lattice together with a “cloud” of magnons. ^{2} Thereby the energy scale of the “coherent” hole propagation is strongly renormalized and is given by the AF superexchange constant $J$. In this way a quasiparticle is formed which is frequently called in the literature a *spin polaron*. ^{3} Actually, such a phenomenon is not only of theoretical importance but also describes a realistic situation found in the photoemission spectra of the parent compounds of the high- ${T}_{c}$ cuprates, as e.g., in ${\text{Sr}}_{2}{\text{CuO}}_{2}{\text{Cl}}_{2}$. ^{4}

A more complex situation can occur in the systems with partly filled degenerate orbitals, where a doped hole may not only couple to magnons but may also couple to the crystal-field excitations. ^{5} Then we may for example arrive at a problem shown schematically in Fig. 1, which is studied in this paper and discussed in more detail below, or end up in two somewhat simpler, though still challenging, situations where at least the coupling to the magnons could be neglected. To get a better insight into the *orbital physics* we give a brief overview of these two latter cases in the next two paragraphs.

First, we recall an example of the $ab$ planes of ${\text{LaMnO}}_{3}$ which has ferromagnetic (FM) spin order and alternating orbital (AO) order of ${e}_{g}$ orbitals in the ground state. ^{6} It has been shown ^{7} that although the hole introduced into such a state does not disturb the FM spin order, it couples to the collective excitations of the AO state (orbitons). Here again a quasiparticle is formed which is called this time an *orbital polaron*. ^{8} It should be noted, however, that the orbital polaron has an even smaller bandwidth than the spin polaron, ^{7} as the orbitons are in general much less mobile than the magnons (or even immobile) due to the lack of the SU(2) symmetry in the orbital systems, ^{9} and almost directional Ising-type superexchange.^{6,10} Actually, one can understand the hole motion in this case in terms of the string picture: ^{11} The hardly mobile orbitons cannot repair the string of the misaligned orbitals in the AO state, which arises by the hole motion. Thus, it is the hole which has to return to the original site and cure the defects by retracing its path, unless it propagated due to small off-diagonal hopping in an ${e}_{g}$ system ^{7} and no defects were created on its path (the latter process also contributes to the above mentioned very small bandwidth of the orbital polaron).

Second, only recently an even more extreme situation of the system with orbital order was investigated:^{12,13} The case of a hole doped into the plane with FM spin order accompanied by the ${t}_{2g}$ AO order [which could correspond not only to the hole introduced into the ordered ground state of ${\text{Sr}}_{2}{\text{VO}}_{4}$ with ${t}_{2g}$ orbitals but also, surprisingly, to those of ${\text{K}}_{2}{\text{CuF}}_{4}$ or ${\text{Cs}}_{2}{\text{Ag}}_{4}$ with ${e}_{g}$ active orbitals (see Ref. 13)]. Also here a quasiparticle is formed due to the dressing of a hole by the collective excitations of the ground state with AO order. However, due to the specific ${t}_{2g}$ orbital symmetries the orbitons are not mobile at all, the off-diagonal hopping is prohibited, and the quasiparticle acquires a finite bandwidth only due to the frequently neglected three-site terms. ^{12} Thus, the string picture dominates the character of the ${t}_{2g}$ orbital polarons even more than in the case of systems with ${e}_{g}$ orbital degrees of freedom.

Finally, after this brief overview, we turn to the problem investigated in the present paper: the motion of a single hole introduced into the ground state with coexisting AF order of $S=1$ spins and ${t}_{2g}$ AO order, shown schematically in Fig. 1. This situation corresponds to a single hole doped into the $ab$ plane of ${\text{LaVO}}_{3}$. ^{14} It differs from all of the three (one—purely spin and the other two—purely orbital) cases described above due to the coexisting AF/AO order. Neither spin nor the orbital background is then transparent for the propagating hole, and the hole has to couple *both* to the magnons and orbitons simultaneously (see Sec. II C). To our knowledge, the mechanism and physical consequences of this coupling have never been studied before. ^{15}

Furthermore, the coexistence of the AO and AF order is an extremely rare case in nature as it formally violates ^{16} the Goodenough-Kanamori rules that state complementary spin and orbital order in the ground state, i.e., either FM spin coexisting with AO order, or AF spin coexisting with FO order—thus, it is interesting to investigate how a hole can propagate in a state which does not fulfill this rule and by its nature contains more quantum fluctuations. Therefore, the main goals of this paper are to verify the following questions: (i) *under which conditions*, if at all, *a quasiparticle can be formed* in such an AF phase with AO order? (ii) If a quasiparticle is formed, *what are its properties* and whether the string picture (which is typical for orbital physics) can explain them? Finally (iii) *what is the role of spins and orbitals in the possible formation of the string* in the spin-orbital polarons? By working out the answers to the above questions we would like to predict the main features of the photoemission spectra of the cleaved samples of ${\text{LaVO}}_{3}$ with polarization corresponding to the $ab$ planes.

Before we move on to present the answers to these questions, let us note another motivation for studying the above problem. It concerns the phase diagram of the lightly doped cubic vanadate ${\text{La}}_{1-x}{\text{Sr}}_{x}{\text{VO}}_{3}$. ^{17} It appears that in this strongly correlated compound the $C$-AF (AF $ab$ planes with FM order along the $c$ direction) coexisting with $G$-AO (AO order in all three cubic directions) order Mott insulating phase is not only stable at $x=0$ but also persists to a relatively high value of doping $x=0.178$. ^{17} Furthermore, the $C$-AF phase remains stable up to an even higher value of $x=0.26$, although in this regime the insulating and orbital ordered phase has already disappeared. ^{17} As in the ionic picture $x$ stands for the fraction of holes doped into the $d$ orbitals of vanadium ions (where a nominal valence upon doping changes as ${d}^{2-x}$), it remains a challenge to explain why the ordered and insulating states persist to such high doping level. ^{18} Besides, somewhat similar phase diagrams have been observed in other doped cubic vanadates, ^{19} such as ${\text{Pr}}_{1-x}{\text{Ca}}_{x}{\text{VO}}_{3}$, ${\text{Nd}}_{1-x}{\text{Sr}}_{x}{\text{VO}}_{3}$, or even to some extent in ${\text{Y}}_{1-x}{\text{Ca}}_{x}{\text{VO}}_{3}$ although in all these cases the lattice (Jahn-Teller and ${\text{GdFeO}}_{3}$-like) distortions contribute significantly to their physical properties, similar to the undoped parent $R{\text{VO}}_{3}$ (with $R=\text{Pr},\text{Nd},\text{Y}$) compounds. ^{20} Here to avoid further complications due to additional orbital interactions imposed by the lattice, we would like to concentrate ourselves on the spin-orbital polarons in (almost) undistorted compound ${\text{LaVO}}_{3}$.

Actually, one should in principle be able to understand some of the features observed in the phase diagram of the doped cubic vanadates by comparing them with those of the high- ${T}_{c}$ cuprates or the colossal magnetoresistive manganites. However, such a comparison only further emphasizes the lack of theoretical understanding of the doped cubic vanadates. First, in the high- ${T}_{c}$ cuprates, such as ${\text{La}}_{2-x}{\text{Sr}}_{x}{\text{CuO}}_{4}$, the AF order disappears very quickly with doping $x$, i.e., already for $x\sim 0.02$. ^{21} This is despite the fact that the value of the superexchange constant $J$ is relatively high there, which naturally leads to larger magnetic energy in the cuprates than in the vanadates. This suggests that, among other factors, it is the orbital dynamics which could be responsible for the totally different behavior of these two classes of compounds upon hole doping. Second, a similar conjecture can be drawn from the comparison between the vanadates and the manganites. In the latter ones, e.g., in ${\text{La}}_{1-x}{\text{Sr}}_{x}{\text{MnO}}_{3}$, the AO orbital Mott insulating state is stable up to $x\sim 0.18$, ^{22} i.e., almost to the same doping level as in ${\text{La}}_{1-x}{\text{Sr}}_{x}{\text{VO}}_{3}$. However, ${\text{La}}_{1-x}{\text{Sr}}_{x}{\text{MnO}}_{3}$ has FM planes but ${\text{La}}_{1-x}{\text{Sr}}_{x}{\text{VO}}_{3}$ has AF planes—this again indicates that it is the orbital dynamics which governs the behavior of holes in the doped cubic vanadates.

Hence, we arrive at the following paradox. On the one hand, it is obvious from the theory that any meaningful study of the doped perovskite vanadates should take into account the spin and orbital degrees of freedom on equal footing. On the other hand, we see from the above conjecture from the experiment that these are the orbital degrees which are, to large extent, responsible for the stability of the observed phases of the lightly doped cubic vanadates. Thus, another aim of the paper is to shed some light on this apparent paradox. Certainly, the studies presented below concern only the case of just a single hole doped into the half-filled compound and, furthermore, they investigate only the situation in the two-dimensional (2D) model focusing on the $ab$ planes. Nevertheless, we hope to understand the generic role of the spin and orbital degrees of freedom in the hole-doped cubic vanadates. Besides, we note that complementary studies presented in Ref. 23 revealed the role of the AO and FM order stable along the third (not studied here) direction in the doped ${\text{La}}_{1-x}{\text{Sr}}_{x}{\text{VO}}_{3}$, and explained the differences between the doped ${\text{Y}}_{1-x}{\text{Ca}}_{x}{\text{VO}}_{3}$ and ${\text{La}}_{1-x}{\text{Sr}}_{x}{\text{VO}}_{3}$ but, by the very nature of that one-dimensional model, could not address the problems mentioned above. We also stress that contrary to the problem solved in this paper, the hole moving along this third $c$ cubic direction couples *only* to *one* type of excitations: either orbitons in the lightly doped $C$-AF phase of ${\text{La}}_{1-x}{\text{Sr}}_{x}{\text{VO}}_{3}$, or magnons in the very lightly doped $G$-AF phase of ${\text{Y}}_{1-x}{\text{Ca}}_{x}{\text{VO}}_{3}$. ^{23}

After a thorough discussion of the motivation standing behind the studies presented in the paper let us now concentrate on the model and the theoretical method used in this work. We use the 2D $t\text{-}J$ model for the ${\text{V}}^{3+}$ ions in the ${t}_{2g}^{2}$ configuration with degenerate ${t}_{2g}$ orbitals and the superexchange given in Refs. 14,24 as *the model* which describes the physical properties of a hole in the $ab$ plane of the lightly doped cubic vanadates. Since any $t\text{-}J$ model which contains Ising-type interactions should be supplemented by the three-site effective hopping terms in order to provide a possibility of coherent hole propagation,^{12,13} we derive these terms and add them to the $t\text{-}J$ model (such a model is then called the strong-coupling model). Next, we solve the strong-coupling model in the case of a single hole doped into the ground state of a Mott insulator by reducing it to the polaronlike model and introducing the self-consistent Born approximation (SCBA). ^{3} Although the final solutions (spectral functions) are obtained numerically on a finite mesh of the momentum $\mathbf{k}$ points, the SCBA method is known to be very accurate to study such problems,^{3,25} and is largely independent of the size of the clusters ^{3} on which the calculations are performed.

The paper is organized as follows. We start with introducing the spin-orbital $t\text{-}J$ model in Sec. II A which is supplemented by the crucial three-site terms for the spin-orbital model, see Appendixes A,B for detailed discussion and derivation. We reduce the derived strong-coupling model ( $t\text{-}J$ model with three-site terms) to the spin-orbital polaron problem in Secs. II B,II C,II D and discuss its parameters in Sec. II E. Next, we formulate the solution of the model using the analytic approach within the SCBA method (see Sec. III). In Sec. IV we present the numerical results and discuss the properties of the quasiparticle states. Finally, in Sec. V we analyze the composite interplay of spin and orbital degrees of freedom in the formation of the spin-orbital polaron whereas in Sec. VI we present our conclusions. We adopt a convention throughout the paper [despite the introductory part (Sec. I) and Appendix A] that the red (blue) color in the online figures denotes features associated solely with spins (orbitals).

##### II. MODEL

###### A. Spin-orbital $t\text{-}J$ model with three-site terms

The Hamiltonian of the spin-orbital $t\text{-}J$ model with the three-site terms (called also the strong-coupling spin-orbital Hamiltonian), relevant for the planes of lightly doped cubic vanadates, consists of three terms,

where the last part ${H}_{3s}$ stands for the three-site terms. The first term in Eq. (2.1), which describes the hopping of the electrons in the constrained Hilbert space, i.e., in the space with singly-occupied or doubly-occupied sites, follows in a straightforward way from the unconstrained hopping of electrons residing in the ${t}_{2g}$ orbitals. The latter one has to fulfill the ${t}_{2g}$ orbital symmetries, meaning that along each bond in the cubic direction only two out of three ${t}_{2g}$ orbitals are active.^{14,26} This means that in the $ab$ plane, which is under consideration here, electrons in ${d}_{yz}\equiv a$ $\left({d}_{zx}\equiv b\right)$ orbitals can hop only along the $b$ $\left(a\right)$ direction. Besides, the ${d}_{xy}$ orbital does not contribute to hopping elements as it lies lower in energy and is always occupied by one electron in the half-filled and lightly hole-doped regime. ^{27} Hence, we arrive at the kinetic $\propto t$ part of Hamiltonian (2.1),

###### (2.2)

Here the constrained operators

mean that the hopping is allowed only in the restricted Hilbert space with not more than one $\left\{a,b\right\}$ electron at each site $\mathbf{i}$ ( $\overline{\sigma}$ stands for the spin component opposite to $\sigma $). Besides, since the Hund’s coupling is large $\left({J}_{H}\u2aa2t\right)$ in the cubic vanadates, ^{28} we project the final states resulting from the electron hopping onto the high spin states, which is denoted by the $\mathcal{P}$ operators in Eq. (2.2).

The middle term in Eq. (2.1) describes the spin-orbital superexchange in cubic vanadates ^{14} —here we use the notation of Ref. 24, see Eqs. (6.5)–(6.7). It reads

where the contribution due to different excitations in the multiplet structure of vanadium ions are ^{29}

Here ${\mathbf{S}}_{\mathbf{i}}$ is a spin $S=1$ operator, ${T}_{\mathbf{i}}^{z}=\left({\stackrel{\u0303}{n}}_{\mathbf{i}b}-{\stackrel{\u0303}{n}}_{\mathbf{i}a}\right)/2$ is a pseudospin $T=1/2$ operator, and $J=4{t}^{2}/U$ is the superexchange constant, with $U$ being the intraorbital Coulomb element on a vanadium ${\text{V}}^{2+}$ ion. As in the spin $t\text{-}J$ model, ^{30} we assume that the hopping $t$ is small, i.e., $t\u2aa1U$ whose assumption is well satisfied by the realistic parameters in the cubic perovskites. ^{14} The factors ${r}_{1}=1/\left(1-3\eta \right)$ and ${r}_{3}=1/\left(1+2\eta \right)$, where

account for the Hund’s coupling ${J}_{H}$ and originate from the multiplet structure of the ${d}^{3}$ excited states characterized by the various allowed spin and orbital configurations. ^{14} Let us note that superexchange Hamiltonian (2.5) was derived ^{24} with the assumption that the ${d}_{xy}$ orbital was singly occupied at each Vanadium ion (see discussion above). Besides, in principle Hamiltonian (2.5) was originally derived for the undoped case and should be modified for the doped case by adding the superexchange interactions due to the existence of the ${d}_{\mathbf{i}}^{1}{d}_{\mathbf{j}}^{1}$ and ${d}_{\mathbf{i}}^{1}{d}_{\mathbf{j}}^{2}$ nearest-neighbor configurations. However, in the discussed here small doping regime, and particularly for a single hole, the modification of magnon and orbiton excitations due to these terms should be very small and we will neglect them.

The last term in Eq. (2.1) denotes the three-site terms which, to our knowledge, have not been derived before for the case of the ${d}^{2}$ systems with spin and orbital degrees of freedom. These terms, although frequently neglected, are crucial for a faithful representation of the Hubbard model in the strong-coupling regime, ^{31} and can play an important role in the coherent hole propagation in orbital systems.^{12,13} However, in the present case the derivation of all possible three-site terms is relatively tedious and leads to a complex expression. Fortunately, it was shown in Ref. 13 (for orbital systems) and in Ref. 32 (for spin systems) that the only three-site terms which occurred to be relevant for the lightly doped systems were those which did not contribute to the coupling between hole and orbital or spin excitation but merely contributed to the free hole motion (see Sec. II D). Nevertheless, we discuss the possible role of the neglected three-site terms in Appendix B. Then, the derived in Appendix A three-site terms which could possibly contribute to the free hole motion read

###### (2.10)

Note that these terms are $\propto J$ and hence are of the same order in ${t}^{2}/U$ as superexchange terms (2.5).

In the next four sections we reduce the model Eq. (2.1) to the polaron problem following Ref. 3. This would enable us to solve the model Eq. (2.1) using the SCBA method.

###### B. Undoped case: low energy excitations

It was shown in Ref. 14 that the ground state of the three-dimensional (3D) spin-orbital superexchange model (2.6) is characterized by the classical phase with $C$-AF spin order and $G$-AO orbital order. Thus, we start from the above *undoped* ground state of the Hamiltonian Eq. (2.1) with classical (Néel) order. As usual, this is not the eigenstate of the Hamiltonian, and thus the full description of the system has to take into account the quantum fluctuations around such a classical ground state. Below, we will calculate them, together with the low-energy (magnon and orbiton) excited states, employing the linear spin wave (LSW) and linear orbital wave (LOW) approximation, following Ref. 33.

First, in the classical state we introduce two sublattices $A$ and $B$ such that: (i) all $a$ $\left(b\right)$ orbitals are occupied in the perfect AO state in sublattice $A$ $\left(B\right)$, and (ii) spins pointing “downwards” (“upwards”) are located on sublattice $A$ $\left(B\right)$, as shown schematically in Fig. 1. Next we rotate spins and orbitals (pseudospins) on sublattice $A$ so that all the spins and pseudospins are aligned in the considered $ab$ plane and take positive values now.

Second, we introduce Schwinger bosons $t$ and $f$ for orbitals and spins such that

with the constraints

Third, we transform the Schwinger boson operators into the Holstein-Primakoff bosons $\alpha $ and $\beta $,

where the above constraints are now no longer needed.

Next, we substitute the above transformations into the Hamiltonian ${H}_{J}$ and skip higher-order terms (using LSW and LOW approximation). The latter approximation physically means that the number of bosons ${\alpha}_{\mathbf{i}}$ and ${\beta}_{\mathbf{i}}$, which describe the fluctuations around the ordered state, is small. This results in the effective substitutions

Finally, we introduce Fourier transformation separately for each sublattice ( $N$ is the total number of sites on both sublattices and the lattice constant $a=1$),

we define operators

and perform the standard Bogoliubov transformation to diagonalize the magnon part of the Hamiltonian

where

with ${\nu}_{\mathbf{k}}=\sqrt{1-{\gamma}_{\mathbf{k}}^{2}}$ and the structure factor for the square lattice ${\gamma}_{\mathbf{k}}=\left(\text{cos}\text{\hspace{0.17em}}{k}_{x}+\text{cos}\text{\hspace{0.17em}}{k}_{y}\right)/2$.

Then, after neglecting constant terms, the LSW and LOW Hamiltonian for magnons and orbitons read

where

and

in agreement with Eq. (6.11) of Ref. 24 and Eq. (11) of Ref. 29. Besides the dispersion relation of the magnons is

where $z=4$ is the coordination number and $S=1$ is the spin value. Let us note that in the regime of reasonable values of $\eta \u220a\left[0,0.20\right]$: ${J}_{O}>0$ and ${J}_{S}>0$, which means that the classical ground state has indeed coexisting AF and AO order. Furthermore, at temperature $T=0$ the considered here classical 2D AO and AF ground state is stable with respect to the quantum fluctuations, both in spin and in orbital channel. In fact, the orbital order is undisturbed by local Ising excitations, while the quantum AF ground state is modified and the order parameter is renormalized with magnon excitations. ^{3}

###### C. Hole coupling with magnons and orbitons

We expect that a doped hole does not modify significantly the classical ground state stable for the half-filled case (see above). This could play a role in the lightly doped regime, but in the case of one hole in the whole plane such a modification is negligible and will be neglected below. Instead, the doped hole may modify its neighborhood by its coupling to the excitations of the classical ground state — magnons and orbitons—which renormalize the hole motion. In order to describe it mathematically, we rewrite ${H}_{t}$ (see below) and ${H}_{3s}$ (see next section) using similar transformations as performed in Sec. II B.

First, we rotate spins and pseudospins on sublattice $A$. Next, we express the electron operators in terms of the Schwinger bosons introduced in Sec. II B,

with the same constraints as in Sec. II B. Note that the factor ${\scriptstyle \frac{1}{\sqrt{2}}}$ is not added “ad hoc” but originates from a detailed check of the validity of the above equations: It should always be present in the case of spin $S=1$ because, e.g., when one annihilates two-boson state with the $f$ operator, then a factor $\sqrt{2}$ appears. Due to this factor and the above constraint on the number of bosons the high spin projection operators $\mathcal{P}$ in ${H}_{t}$ are no longer needed (i.e., quantum double exchange ^{34} factors are implicitly included in this formalism).

Next, similarly as in Sec. II B, we transform the Schwinger bosons into the Holstein-Primakoff bosons, skip all terms containing more than two bosons, perform Fourier transformation for bosons and (additionally) for holons here, introduce ${\alpha}_{\mathbf{k}\pm}$ operators, and finally perform Bogoliubov transformation to arrive at the Hamiltonian

###### (2.39)

where

with $\mu =x,y$ and $\overline{\mathbf{k}}=\mathbf{k}-{\mathbf{q}}_{1}-{\mathbf{q}}_{2}$; the coefficients $\left\{{u}_{{\mathbf{q}}_{1}},{v}_{{\mathbf{q}}_{1}}\right\}$ are standard Bogoliubov factors (2.30).

###### D. Doped hole: The free dispersion

After performing similar transformations as the ones introduced in Sec. II C we obtain that the three-site terms, Eq. (2.10), lead to the following Hamiltonian for the holes:

where

and the free dispersion relations for two orbital flavors are

Note that we neglected all of the three-site terms which lead to the coupling between holes and magnons and/or orbitons. This approximation is physically well justified since all such terms would be of the order of $J/4$, i.e., much smaller than the terms in Eq. (2.39). See also Appendix B for further discussion.

###### E. Polaron model and its parameters

The above considerations demonstrate that in the lightly doped case, when the classical spin and orbital ordered ground state present in the half-filled case is robust, strong-coupling model (2.1) can be reduced to an effective model

[see Eqs. (2.31) and (2.39),(2.40),(2.41)]. Actually, this is a polaron-type model with the coupling between fermions (holes) and bosonic excitations (orbitons and magnons), which is relatively easy to solve (see below). The validity of the mapping between the two models was thoroughly discussed in Refs. 3,13.

Note that original strong-coupling model (2.1) has three parameters $\left\{J,\eta ,t\right\}$, whereas the effective polaron model given by Eq. (2.45) is more conveniently analyzed using instead four parameters $\left\{{J}_{S},{J}_{O},\tau ,t\right\}$, which determine the scale of spin and orbital excitations, as well as free hole propagation (due to the three-site terms) and the vertex function $\left(t\right)$ (see below). In what follows we will use either one of these two parameter sets (or even both of them), depending on the context. Hence, we plot in Fig. 2 the functional relation between the parameters $\left\{{J}_{S},{J}_{O},\tau \right\}$ on Hund’s exchange $\eta $. While the magnon energies $\propto {J}_{S}$ decrease with $\eta $, the energy scale of orbitons $\propto {J}_{O}$ increases rapidly, so the latter excitations are expected to play an important role in the realistic regime of parameters.

##### III. GREEN’S FUNCTIONS AND SPECTRAL FUNCTIONS

###### A. Self-consistent Born approximation

The spectral properties of the hole doped into the AF/AO ground state $|{\Phi}_{0}\u27e9$ with energy ${E}_{0}$ of the strong-coupling model Eq. (2.1) at half-filling follow from the Green’s functions,

However, due to the mapping of the strong-coupling model onto the polaron model performed in Secs. II B,II C,II D,II E, it is now convenient to express the above Green’s functions in terms of the polaron Hamiltonian $\mathcal{H}$. This requires that we first write down the electron operators in terms of the operators used in Eq. (2.45),

Second, the ground state $|{\Phi}_{0}\u27e9$ is now a physical vacuum $|0\u27e9$, with respect to the Bogoliubov operators ${\stackrel{\u0303}{\alpha}}_{\mathbf{k}\pm}$ and the operators ${\beta}_{\mathbf{k}}$, with energy $E$. Then, we arrive at the following relations:

whereas

and we skip the two latter Green’s functions in what follows. Note that the above set of equations follows from the fact that ${\beta}_{\mathbf{k}}|0\u27e9=0$ and the inequalities are due to

where ${n}_{0}$ is the average number of spin deviations in the 2D ground state $|0\u27e9$. Note that below we will eliminate the ground-state energy $E$ to simplify equations.

As seen in Sec. II C the vertices in the spin-orbital model are slightly more complex than for the standard spin $t\text{-}J$ model: ^{3} (i) one always finds two boson and two holon lines at each vertex (instead of just one boson and two holon lines), (ii) the two sublattice structure plays an essential role (this resembles the orbital case), and (iii) one has two kinds of magnons (which we had to introduce in order to keep track of the lattice index). Also, in the LSW and LOW order all the terms $\propto {\alpha}_{\mathbf{i}}^{\u2020}{\beta}_{\mathbf{i}}$ do not contribute (the self-energies for them would require altogether four boson lines instead of just two for $\propto {\alpha}_{\mathbf{i}}{\beta}_{\mathbf{i}}$). ^{35} Hence, we largely follow the route proposed for the orbital ${t}_{2g}$ model ^{13} and obtain the following SCBA equations for the self-energies (see also Fig. 3):

###### (3.14)

###### (3.15)

where the factor 2 in front of each vertex comes from the fact that we have two kinds of magnons [ $+$ and $-$, see Eq. (2.28)] and hence two distinct diagrams. Fortunately, this factor cancels with one of the 2’s in the denominator and we obtain that the coupling constant is ${\left(t/\sqrt{2}\right)}^{2}$, i.e., we recover this factor $1/\sqrt{2}$ which originates from the quantum double exchange. ^{34} The above equations should always be supplemented by the Dyson’s equations: (Fig. 3),

which together form a self-consistent set of equations which can be solved numerically by iteration.

Finally, we can calculate the spectral functions for a hole created in $a$ and $b$ orbital,

###### (3.18)

###### (3.19)

where we introduced a factor of 2 in front of the definition of the spectral function ${A}_{\gamma}\left(\mathbf{k},\omega \right)$ for convenience.

###### B. Analytic structure and quasiparticle pole

It occurs that in the case when the three-site terms are absent (i.e., at $\tau \equiv 0$) one can easily prove two important properties of the SCBA Eqs. (3.14) and (3.15): (i) the self-energies are $\mathbf{k}$ independent and (ii) the spectral functions contain the quasiparticle state for finite value of the exchange parameter $J$ (equivalently ${J}_{S}$ *or* ${J}_{O}$).

First, we show property (i). Since we assumed that $\tau =0$, we can rewrite SCBA Eqs. (3.14) and (3.15) together with Dyson’s Eqs. (3.16) and (3.17) in the following manner:

###### (3.20)

###### (3.21)

which after substitution ${\mathbf{q}}_{2}\to \mathbf{k}-{\mathbf{q}}_{2}$ in the sums leads to

###### (3.22)

###### (3.23)

where we have used a short-hand notation

Since the self-energies on the right-hand side of Eqs. (3.22) and (3.23) do not depend on $\mathbf{k}$ we are allowed to drop the momentum dependence of the self-energies. Hence, the spectral functions are also momentum independent in the case of $\tau =0$, and we recognize that the dependence of the spectral functions on $\mathbf{k}$ may originate only from the three-site terms.

Second, using the dominant pole approximation ^{2} we show that the quasiparticle state exists [property (ii)] if $J$ is finite (i.e., ${J}_{S}$ *or* ${J}_{O}$ are finite). Hence, following Kane *et al.*, ^{2} we assume that the Green’s function can be separated into the part containing the pole and the part responsible for the incoherent processes:

where

and the pole positions:

have to be determined self-consistently following Eqs. (3.26) and (3.27).

Next, following Ref. 2, it is straightforward to derive the upper bound for the residues $\left\{{a}_{A},{a}_{B}\right\}$:

If the sums in the above equations are divergent, then the upper bounds for the residues are equal zero and the Green’s functions do not have the quasiparticle pole. Hence, we need to check the behavior for small values of the momenta ${\mathbf{q}}_{1}$. Then ${\omega}_{{\mathbf{q}}_{1}}\propto \left|{\mathbf{q}}_{1}\right|$, but e.g., ${\left({u}_{{\mathbf{q}}_{1}}{\gamma}_{{q}_{2y}-{q}_{1y}}+{v}_{{\mathbf{q}}_{1}}{\gamma}_{{q}_{2y}}\right)}^{2}\sim \left|{\mathbf{q}}_{1}\right|{\left({\gamma}_{{q}_{2y}}-{\widehat{q}}_{1}\cdot \nabla {\gamma}_{{q}_{2y}}\right)}^{2}$. Thus, if at least either ${J}_{S}$ or ${J}_{O}$ is finite, then there are no divergences in Eqs. (3.32) and (3.33). Consequently, under the same conditions the quasiparticle state exists.

##### IV. NUMERICAL RESULTS

###### A. Spectral functions

We calculated the spectral functions ${A}_{a}\left(\mathbf{k},\omega \right)$ and ${A}_{b}\left(\mathbf{k},\omega \right)$ by solving SCBA Eqs. (3.14) and (3.15) on a mesh of $16\times 16$ $\mathbf{k}$-points. The results for two different values of the superexchange constant $J=0.2t$ and $J=0.6t$ are shown in Figs. 4 and 5, respectively. Also we assume the Hund’s coupling to be close to a realistic value in ${\text{LaVO}}_{3}$, i.e., $\eta =0.15$. ^{36} We emphasize that such a value of the Hund’s coupling is not only realistic but together with the observed value of $J=0.2t$ gives a reasonable value of the spin-only exchange constant ${J}_{S}=0.03t$ (which is in agreement with the observed value of the Néel temperature in ${\text{LaVO}}_{3}$ ${T}_{N}\simeq 143\text{\hspace{0.5em}}\text{K}$),^{17,37} cf. Sec. VIB of Ref. 24 for a more detailed discussion on this issue. Thus, the spectral functions shown in Fig. 4 are calculated for the realistic values of parameters in ${\text{LaVO}}_{3}$ whereas those shown in Fig. 5 will serve as a comparison with the regime of large $J$.

Let us first discuss the results presented in Fig. 4. The quasiparticle peak in the low energy part of the spectrum is clearly visible and confirms the analytic predictions presented in Sec. III B. However, the quasiparticle state has a rather weak one-dimensional (1D) dispersion: along the ${k}_{x}$ direction for holes doped into the $b$ orbitals and along the ${k}_{y}$ direction for holes doped into the $a$ orbital. Since such a dispersion was not present when the three-site terms were neglected (see Sec. III B), we can immediately ascribe the onset of the 1D dispersion in the spectra to the three-site terms. Indeed, these terms are responsible for the coherent hole propagation in the ground state with AO order and the hole dispersion is dressed by incoherent processes. This phenomenon is quite well understood in the framework of the orbital-only strong-coupling model—we refer to Refs. 12,13 for the more detailed discussion.

Besides, the excited states form a ladderlike spectrum, also with a small 1D dispersion. Similarly to the quasiparticle state the excited states resemble qualitatively the spectral functions calculated for the purely orbital model (see Refs. 12,13). However, here we see that the spin-orbital spectra are quantitatively different than the orbital ones. In particular, they are quite similar to those obtained for $J=0.4t$ in the purely orbital model (see Fig. 9 in Ref. 13) although this is just a coincidence since all of the exchange constants used in the spin-orbital model ( $J$, ${J}_{S}$, or ${J}_{O}$) are much smaller than $0.4t$.

Looking at the spectral functions shown in Fig. 5 we note that the ladderlike spectrum almost disappears and the quasiparticle dispersion increases to $\sim 0.4t$ when the value of the superexchange constant is increased to $J=0.6t$. However, the qualitative features in the dispersion of the quasiparticle peak stay mostly unchanged and the 1D character of the dispersion relation is preserved. Quantitatively, both the quasiparticle spectral weight and the bandwidth increase quite drastically.

To conclude, the spin-orbital spectral functions form ladderlike spectra with a small 1D dispersion and have many similarities with the purely orbital spectra of the ${t}_{2g}$ model. ^{12} Still, however, there are few relatively important differences with the orbital model which suggest that the spin-orbital case is more complex and the quasiparticle behavior more subtle than for the purely orbital model. For the understanding of this problem we refer to Sec. V whereas in Sec. IV B we further investigate the properties of the quasiparticle arising in the spectral functions of the spin-orbital model.

###### B. Quasiparticle properties

In this section we analyze in more detail the quasiparticle properties of the spectral functions of the spin-orbital model, cf. Fig. 6. In addition, we compare these properties with those calculated for the spin strong-coupling model of Ref. 32 and the orbital strong-coupling model of Ref. 12.

By looking at the results for the spin-orbital model in Fig. 6 we see that: (i) the quasiparticle bandwidth $W$ is strongly renormalized from its respective free value (which is $W=8\tau $ as it may only originate from the three-site term dispersion relation, see Secs. III B and IV A) and is proportional to ${J}^{2}$ for small $J$ $\left(J\lesssim 0.6t\right)$ and to $J$ in the regime of large $J$ $\left(J\gtrsim 0.6t\right)$; (ii) the quasiparticle spectral weight ${a}_{\text{QP}}$ grows considerably with increasing $J$; and (iii) the pseudogap $\Delta $ (the energy distance between the quasiparticle state and the first-excited state) exists and roughly scales as $t{\left(J/t\right)}^{2/3}$ although for larger values of superexchange $J$ some deviations from this law were observed.

Thus, on the one hand, we see that the quasiparticle properties of the spin-orbital model are qualitatively distinct than those of the spin model. We should stress that although the spin model used here is the $S=1/2$ $t\text{-}J$ model with the three-site terms ^{32} but we expect similar generic behavior in the case of spin $S=1$ (see also Sec. V A). Hence, as there is no qualitative agreement with the model for spin $S=1/2$ there would neither be a qualitative agreement with the model for spin $S=1$.

On the other hand, we note remarkable similarities between the present model and the orbital model (see also Fig. 6); the generic behavior of the bandwidth, the spectral weight, and of the pseudogap seems to be almost identical in both models. There are, however, two important differences. First, the $t{\left(J/t\right)}^{2/3}$ law does not describe the behavior of the spin-orbital pseudogap as well as in the orbital case. Second, if we assume that the spin-orbital model is just qualitatively similar to the orbital model, then it should be possible to rescale all the quasiparticle properties with such an effective value of the superexchange $J$ that they coincide with the results for the orbital model. This, however, is not possible. For example, from the fits to the $t{\left(J/t\right)}^{2/3}$ law we can deduce that such an effective value of the superexchange would be ${J}_{\text{eff}}={\left({a}_{SO}/{a}_{O}\right)}^{2/3}=0.79J$, where ${a}_{SO}=1.94$ $\left({a}_{O}=1.66\right)$ is the fitted coefficient which multiplies the $t{\left(J/t\right)}^{2/3}$ law in the spin-orbital (orbital) case. This single parameter does not suffice, however, as a similar fit to the quasiparticle spectral weight would require such a ${J}_{\text{eff}}$ to be bigger than $J$.

In conclusion, the properties of the quasiparticle state in the spin-orbital model resemble those found for the quasiparticle in the purely orbital ${t}_{2g}$ model.^{12,13} However, a more detailed analysis presented in this section reveals that such a correspondence is rather “superficial” and that few of the qualitative features of both models are different. We refer to Sec. V for the thorough understanding of this phenomenon.

##### V. ORIGIN OF SPIN-ORBITAL POLARONS

###### A. Comparison with spin and orbital polarons

The main task of this and the following two sections is to understand the spectral functions and the quasiparticle properties of the spin-orbital model. In particular, we not only want to understand the peculiar similarities or rather small differences between the spin-orbital model and its purely orbital counterpart, which were discussed in the last two sections (see Secs. IV A and IV B), but also we want to understand why the spin physics seems to be “hidden” in the present spin-orbital system. Note that for the sake of simplicity in this section we will entirely neglect the three-site terms [Eq. (2.10)] in spin-orbital model (2.1) (or in other words we will put $\tau \equiv 0$). This is motivated by the fact that the role of the three-site terms is solely to provide dispersion to the spectra (see Secs. III B and IV B) and the origin and character of these dispersive features is quite well understood by now (see. Ref. 13).

We start the analysis of spectral functions by introducing two artificial toy models in order to extract the features related to spin and orbital excitations separately, whose results will be later compared to those of the $t\text{-}J$ spin-orbital model. First, we define the following $t\text{-}{J}_{S}$ toy-spin model:

###### (5.1)

where spin operators $\left\{{\mathbf{S}}_{\mathbf{i}}\right\}$ stand for $S=1$ spins, the constrained operators, ${\stackrel{\u0303}{c}}_{\mathbf{i}\sigma}^{\u2020}={c}_{\mathbf{i}\sigma}^{\u2020}\left(1-{n}_{\mathbf{i}\overline{\sigma}}\right)$, exclude double occupancies from the Hilbert space, and the operators $\mathcal{P}$ project onto the high spin states. Note that here the superexchange energy scale is ${J}_{S}$ [see Eq. (2.35)] and not $J$. Hence, it is defined in such a way that it mimics the influence of AO order on the spin subsystem. On the other hand, the kinetic energy is blind to the AO here and Eq. (2.2) reduces to the kinetic part of Eq. (5.1) only if the orbitals form an orbital liquid state. This is an obvious logical inconsistency but our aim here is to see what happens when the joint spin-orbital dynamics in the kinetic energy is entirely neglected. It is also the reason why we call model (5.1) the toy-spin model.

Second, we define the following $t\text{-}{J}_{O}$ toy-orbital model:

###### (5.2)

where pseudospin $T=1/2$, and the constrained operators ${\stackrel{\u0303}{b}}_{\mathbf{i}}^{\u2020}={b}_{\mathbf{i}}^{\u2020}\left(1-{n}_{\mathbf{i}a}\right)$ and ${\stackrel{\u0303}{a}}_{\mathbf{i}}^{\u2020}={a}_{\mathbf{i}}^{\u2020}\left(1-{n}_{\mathbf{i}b}\right)$ exclude double occupancies. Similarly as in the toy-spin model defined above, the superexchange energy scale is not $J$ but ${J}_{O}$ [see Eq. (2.34)], which mimics the influence of spin system with the AF order on the orbital one. Also, the kinetic energy is blind now to the AF order and Eq. (2.2) reduces to the kinetic part of Eq. (5.2) only if the spins form the FM order. This again is logically inconsistent, but we wish to consider this situation for the same reasons as for the spins (see above).

Next, we solve these models using the SCBA method in the case of the one hole doped into the AF (AO) state for the toy-spin (toy-orbital) model. We do not show here the respective SCBA equations as these are straightforward and follow from those written in Ref. 3 (Ref. 13) in the case of the spin (orbital) model. One only has to remember to substitute in the respective SCBA equations $J\to -{J}_{S}$, $S\to 1$, and (due to the quantum double exchange factor) also $t\to t/\sqrt{2}$ in Ref. 3 in the spin case, and $J\to -{J}_{O}$, ${E}_{0}\to 0$, and $\tau \to 0$ in Ref. 13 in the orbital case. We then calculate the respective spectral functions on a mesh of $16\times 16$ $\mathbf{k}$ points. Note that since $\tau \equiv 0$ the spectral functions for both orbital flavors are equal, i.e., ${A}_{a}\left(\mathbf{k},\omega \right)={A}_{b}\left(\mathbf{k},\omega \right)\equiv A\left(\mathbf{k},\omega \right)$.

Finally, we compare the results obtained for the above toy models with those obtained for the spin-orbital $t\text{-}J$ model, Eq. (2.2) and (2.5). We show the results for two different values of $J=0.2t$ and $0.6t$ (see Figs. 7 and 8). In addition, we calculate the results for two different values of the Hund’s coupling $\eta =0$ and $\eta =0.15$ (see left and right panels of each figure).

Let us first look at the physical regime of finite $\eta =0.15$ and $J=0.2t$ [see Figs. 7(d)–7(f)]. On the one hand, the spin-orbital spectral function [panel (d)] resembles qualitatively the ladder spectrum found in the orbital model [panel (f)] although the quantitative comparison reveals strong differences between the two models. On the other hand, the spin-orbital spectral function is entirely different from the $\mathbf{k}$-dependent spin spectral function [panel (e)]. One also finds a somewhat similar behavior for the case of $\eta =0.15$ and $J=0.6t$ [see Figs. 8(d)–8(f)]. Here, however, the spin-orbital spectrum is much more different than the orbital spectrum so the renormalization by the spin part is stronger.

Even more inquiring behavior is found in the nonphysical regime of $\eta =0$ (which, however, is an interesting mathematical limit). ^{36} Then neither of the panels shown in Figs. 7(a)–7(c) or Figs. 8(a)–8(c) is similar to each other. This means that even the orbital model is entirely different in this regime than the spin-orbital model. This is because in this limit the hole moves in the orbital model incoherently as ${J}_{O}=0$ for $\eta =0$ [see e.g., Fig. 7(c)]. However, apparently in the spin-orbital model with ${J}_{O}=0$ and small but finite ${J}_{S}$ the hole moves in stringlike potential, see e.g., Fig. 7(a). This means that the onset of the ladderlike spectrum in the spin-orbital model in this regime cannot be explained easily in terms of the purely orbital model.

In addition, one sees that whereas for different values of $\eta $ but the same values of $J$ one gets rather similar spin-orbital spectra, the spectra found for the toy-orbital model are different. On the contrary, for different values of $J$ and the same values of $\eta $ the spin-orbital spectra are rather distinct but the spin spectra do not change much. This is another argument which suggests that neither the toy-spin model nor the toy-orbital model can explain alone the features found in the spin-orbital spectra.

To conclude, we note that the joint spin-orbital dynamics in the kinetic part of the spin-orbital model plays a significant role for the coherent hole motion. When spin and orbital degrees of freedom are separated from each other (disentangled), the purely spin or orbital models cannot reproduce the spectral function found for the spin-orbital $t\text{-}J$ model. Hence, indeed the spin-orbital spectral functions resemble the orbital ones only superficially and it is the peculiar interplay of the spins and orbitals, studied in Sec. V B, which leads to the calculated spectra.

###### B. Spin-only string excitations at $\eta =0$

In this section we attempt ^{36} to understand the spin-orbital spectra in the limit of $\eta =0$ by assuming that the spins $S=1$ are purely classical and interact by Ising superexchange. Hence, we skip the spin-flip terms $\propto {S}_{\mathbf{i}}^{+}{S}_{\mathbf{j}}^{-}$ in Hamiltonian (2.5), rewrite SCBA Eqs. (3.14) and (3.15) in this case, and finally try to compare the spectral functions calculated in this regime with the ones obtained for the full model, given by Eq. (2.1). In addition, we not only assume $\eta =0$ (which implies ${J}_{O}=0$ in particular) but also $\tau \equiv 0$ (similarly as in the previous section). In Sec. V C we discuss the impact of the finite value of these parameters on the results obtained under these assumptions.

Since ${u}_{{\mathbf{q}}_{1}}=1$, ${v}_{{\mathbf{q}}_{1}}=0$, and ${\omega}_{{\mathbf{q}}_{1}}=zS$ for the Ising superexchange, ^{3} we can rewrite self-energy Eqs. (3.20) and (3.21),

where we already substituted ${\mathbf{q}}_{1}\to \mathbf{k}-{\mathbf{q}}_{1}-{\mathbf{q}}_{2}$ in the sums. Then the above self-consistent equations are momentum independent and we obtain

since $1/N{\sum}_{{\mathbf{q}}_{1}}{\gamma}_{{q}_{1\nu}}^{2}=1/\left(2z\right)$, where $\nu =x,y$.

We have solved Eqs. (5.5) and (5.6) self-consistently. The respective spectral functions are shown in Fig. 9 for $J=0.2t$ (i.e., ${J}_{S}=0.05t$) and $J=0.6t$ (i.e., ${J}_{S}=0.15t$). As expected, one obtains a typical ladder spectrum in the entire regime of $J$. However, we see that the results resemble those obtained for full spin-orbital model (2.1) with $\eta =0$ and $\tau =0$: Although the spectrum of the full spin-orbital model contains some incoherent part, the ladder peaks of the full spin-orbital model and of its classical version (with Ising spin superexchange) almost coincide. In addition, the incoherent bandwidth in the $J=0$ limit is $W=4t$ in the classical case whereas it is only slightly reduced in the quantum model $\left(W\sim 3.7t\right)$. ^{38} For finite $J$ this results in the small shift of the peaks in the full spin-orbital model with respect to its classical counterpart. Altogether, this suggests that the classical and the full (quantum) versions of the spin-orbital models are to a large extent equivalent.

We note that Eqs. (5.5) and (5.6) are almost identical to the SCBA equations for the hole moving in the $S=1/2$ spin Ising model [cf. Eq. (20) in Ref. 3]. The only differences are: the self-consistent dependence of the self-energies on two different sublattices, the reduction in the nearest neighbors by a factor 1/2 (in the numerator), the reduction factor $1/\sqrt{2}$ for the hopping, and the increase in the magnon excitation energy by a factor of 2. Whereas the first two imply the zigzag hole motion in the ordered state, ^{13} the two others merely mean that the hole moves in the spin $S=1$ system, cf. similar factors for the triplet hole in the $t\text{-}J$ model. ^{39} This demonstrates that the hole motion in the full spin-orbital model with $\eta =0$ and $\tau =0$ can be quite well approximated by the zigzag hole motion in the $S=1$ spin Ising model.

The whole analysis presented above leads us to the conclusion that in the limit of ${J}_{0}=0$ and $\tau =0$ the hole moves in the spin and orbitally ordered plane in the following way: (i) the orbitals force the hole to move along the zigzag paths even in the limit of ${J}_{0}=0$, (ii) the orbitals force the hole to retrace its path again even in the limit of ${J}_{0}\to 0$—this resembles the situation discussed by Brinkman and Rice (Ref. 40) for the spin system, where the hole in the Ising model with $J=0$ always has to retrace its path, (iii) the coherent hole motion by the coupling to the spin fluctuations cannot occur in the ground state as then the hole would not retrace its path, so (iv) instead the *hole creates strings in the spin sector* which are erased when the hole retraces its path. Thus, we note here a complex interplay of spins and orbitals. In particular, due to point (ii) *the orbitals force the spins to effectively act on the hole as classical object*.

###### C. Composite string excitations at $\eta >0$

One may wonder how to extend the above understanding of the role of spin and orbital excitations in the formation of spin-orbital polarons to the case of finite values of orbital exchange interaction ${J}_{O}>0$ at $\eta >0$ (see Fig. 2), or in the presence of three-site terms at finite $\tau $. On the one hand, including the nonzero value of ${J}_{O}$ merely leads to the substitution of ${J}_{S}zS\to {J}_{S}zS+{J}_{O}$ in Eqs. (5.3) and (5.4) and consequently also in Eqs. (5.5) and (5.6); this increases the stringlike potential acting on the hole, coming in this case also from the orbital sector. On the other hand, including the three-site term $\tau >0$ results in the energy shift $\omega \to \omega +{\epsilon}_{A}\left(\mathbf{k}\right)$ and $\omega \to \omega +{\epsilon}_{B}\left(\mathbf{k}\right)$ in Eqs. (5.3) and (5.4) which results in the $\mathbf{k}$ dependence of the spectral functions [see Eqs. (5.5) and (5.6)]. However, we can still solve the model using the SCBA. The results (not shown) resemble those found in Figs. 4 and 5: it is again only the incoherent part which is slightly enhanced in full spin-orbital model (2.1) whereas in its classical counterpart it is suppressed. Furthermore, all of the quasiparticle properties of the full spin-orbital model shown in Fig. 6 are *almost perfectly* reproduced by the classical spin-orbital model (not shown)—with the only slight discrepancy occurring in the region of the slight deviation from the $t{\left(J/t\right)}^{2/3}$ law in the pseudogap of the full spin-orbital model.

Finally, there remains just one subtle issue: how to explain the appearance of the small momentum-independent incoherent part in the high-energy part of the spectrum for the quantum spin-orbital model, cf. Fig. 9 in the case of $\eta =0$ [a similar incoherent dome is visible in the spectrum for finite $\eta $ (not shown)]. This can be understood in the following way: although the hole has to return to the original site (due to the orbitals) the magnons present in the excited states can travel freely in the system. Hence, the energies of the excited states can no longer be classified merely by the length of the retraceable paths (as it would be the case in the classical model with no dispersion in the magnon spectrum). This results in the small incoherent spectrum which surrounds each of the peak of the ladder spectrum visible in for example Fig. 9. Furthermore, this incoherent spectrum grows with increasing ${J}_{S}$ as then the velocity of the magnons increases.

To conclude, we note that the quantum spin fluctuations are to a large extent suppressed in the spin-orbital model by the simultaneous coupling of the hole to *both* spin and orbital excitations. In particular, they do not affect the quasiparticle state and merely add as a small incoherent spectrum in the incoherent high-energy part of the ladder spectrum. This is due to the classical character of the orbitals which confine the hole motion and prohibit its coherent motion by the coupling to the quantum spin fluctuations. On the other hand, the hole still couples to the spin degrees of freedom, mostly, in a classical way, i.e., by generating string potential due to defects created by hole motion. Thus, *the string* which acts on the hole moving in the plane with AO and AF order *is of the composite orbital and spin character*. This not only explains the peculiar correspondence between the orbital and spin-orbital model but also explains that the spins “do not hide behind the orbitals” but play an active role in the lightly doped spin-orbital system.

##### VI. CONCLUSIONS

Let us now summarize the results of the present study. We studied the motion of single hole doped into plane of the Mott insulator ${\text{LaVO}}_{3}$, with coexisting AF and AO order, ^{14} shown in Fig. 1, which violates the Goodenough-Kanamori rules. The framework to describe the hole motion in this situation is the respective $t\text{-}J$ model with spin-orbital superexchange,^{14,24} supplemented by the three-site terms required here ^{12} (which were derived in Appendixes A,B), and we used this model as the starting point of our analysis. Next, we reduced this model to the polaron-type model and showed that the hole couples in this case simultaneously to the collective excitations of both the AF state (magnons) and the AO state (orbitons). The problem of finding the spectral function was formulated and solved within the SCBA method. It should be emphasized that, to the best of our knowledge, the problem of a hole coupled *simultaneously* to the magnons and orbitons has never been addressed before. ^{15}

We would like to note that the presented spectral functions, such as the ones shown in Fig. 4, predict the main features to be measured in future photoemission experiments on the cleaved ${\text{LaVO}}_{3}$ sample. The reader may wonder whether (apart from the matrix-elements effects responsible for certain redistribution of spectral intensity) any other processes, such as for example the electron-phonon interaction, would affect hole motion to such an extent that the spectral functions calculated here would change qualitatively. Although we have not made any calculations for such a more complex case so far, we suggest that they would only enhance the ladder spectrum obtained here, since typically the studied mechanisms are only responsible for further localization of the hole. ^{8}

Next, it was calculated both analytically using the dominant pole approximation and numerically using the SCBA equations that the spectral functions contain *a stable quasiparticle peak, provided the value of the superexchange * $J$ was finite. We emphasize that this result was nontrivial as it was not clear *a priori* whether the coupling between a hole and two excitations would lead to a stable quasiparticle state—e.g., the coupling between the hole and two magnons does not lead to the stable quasiparticle peak, cf. Ref. 32 and Appendix B. However, since the orbitons are massive excitations, the hole does not scatter strongly on them and the quasiparticle solution survives.

Furthermore, we studied in detail the properties of the quasiparticle states. In particular, we looked at the differences between the well-known spin ^{3} or orbital^{7,12,13} polarons and the obtained here spin-orbital polarons. We checked that all of the typical quasiparticle properties of the spin-orbital polarons such as the bandwidth, the quasiparticle spectral weight, and the pseudogap (the distance between the quasiparticle peak and the next excited state) are qualitatively similar to those of the ${t}_{2g}$ orbital polarons studied recently,^{12,13} and *it is the string picture which dominates in the quasiparticle properties*. For example, the bandwidth scales as $\propto {J}^{2}$ as a result of the renormalization of the three-site terms $\propto J$, similarly as in the purely orbital ${t}_{2g}$ model.^{12,13} We suggest that the occurrence of the small but still finite bandwidth confirms the idea of the absence of hole confinement in transition-metal oxides with orbital degeneracy presented in Ref. 12.

Finally, a detailed investigation leads us to the conclusion that the spin-orbital polarons are microscopically much more complex than either spin or orbital polarons, and resemble the orbital polarons only superficially. Actually, we have shown that the spin degrees of freedom also play a significant role in the formation of the spin-orbital polarons, although their dynamics is quenched and they are forced by the orbitals to act on the hole as the classical Ising spins. This happens because the orbitals confine the hole motion by forcing the hole to retrace its path which means that the hole motion by its coupling to the quantum spin fluctuations is prohibited. Thus, *the stringlike potential which acts on the hole is induced by the orbitals although it has a joint spin-orbital character*. Furthermore, it occurred that it is only in the excited part of the spectrum that the quantum spin fluctuations contribute and are responsible for a very small incoherent dome in the spectral function.

There is also another important feature of this orbitally induced string formation: it could be understood as a topological effect. This is due to the fact that it happens even if the energy of the orbital excitations is turned to zero, i.e., when the hole moves in the orbital sector incoherently. Hence, the mere presence of orbitals is sufficient to obtain the (almost) classical ladder spectrum of a hole doped into the AF and AO ordered state. This result, in connection with the fact that the mother-compound of the superconducting iron-pnictides shows a variety of spin-orbital phenomena, ^{41} suggests that further investigation of the hole propagation in spin-orbital systems is a fascinating subject for future studies.

Having explained all the main goals of the paper, posed in the sixth paragraph of Sec. I, we can now try to explain the paradox mentioned in the second part of the Introduction. Actually, in the last two paragraphs above we noted that the spin-orbital polarons resemble superficially the orbital polarons but their nature is of joint spin and orbital character. On the condition that this result could be also extended to the case of the finite doping we immediately solve the paradox that experimentally the orbitals play a significant role in the doped cubic vanadates, although theoretically we should take into account both spin and orbital degrees of freedom. Certainly, any straightforward extension of the present one-hole result to the regime of finite doping has to be rather speculative, so further theoretical studies on the doped cubic vanadates are needed to verify our explanation of this paradox.

## References

- L. N. Bulaevski, E. L. Nagaev, and D. I. Khomskii, Sov. Phys. JETP
**27**, 836 (1968). - C. L. Kane, P. A. Lee, and N. Read, Phys. Rev. B
**39**, 6880 (1989). - G. Martínez and P. Horsch, Phys. Rev. B
**44**, 317 (1991). - A. Nazarenko, K. J. E. Vos, S. Haas, E. Dagotto, and R. J. Gooding, Phys. Rev. B
**51**, 8676 (1995); A. Damascelli, Z. Hussain, and Z.-X. Shen, Rev. Mod. Phys.**75**, 473 (2003). - J. Zaanen and A. M. Oleś, Phys. Rev. B
**48**, 7197 (1993). - L. F. Feiner and A. M. Oleś, Phys. Rev. B
**59**, 3295 (1999). - J. van den Brink, P. Horsch, and A. M. Oleś, Phys. Rev. Lett.
**85**, 5174 (2000). - Apart from considered here processes connected with the hopping $t$ (string formation) other coupling mechanisms, e.g., due to the electron-phonon interaction, may contribute to orbital polaron formation; see, e.g., R. Kilian and G. Khaliullin, Phys. Rev. B
**60**, 13458 (1999). - L. F. Feiner and A. M. Oleś, Phys. Rev. B
**71**, 144422 (2005). - J. van den Brink, P. Horsch, F. Mack, and A. M. Oleś, Phys. Rev. B
**59**, 6795 (1999). - Note that although the string picture alone cannot explain the previously mentioned hole motion in the AF ordered state, it may serve as a perfect starting point for the investigation of the behavior of holes doped into the AF phases; see P. Wróbel, W. Suleja, and R. Eder, Phys. Rev. B
**78**, 064501 (2008). - M. Daghofer, K. Wohlfeld, A. M. Oleś, E. Arrigoni, and P. Horsch, Phys. Rev. Lett.
**100**, 066403 (2008). - K. Wohlfeld, M. Daghofer, A. M. Oleś, and P. Horsch, Phys. Rev. B
**78**, 214423 (2008). - G. Khaliullin, P. Horsch, and A. M. Oleś, Phys. Rev. Lett.
**86**, 3879 (2001). - A simpler case of a hole doped into the generic spin-orbital $t\text{\u2212}J$ model was studied in Ref. 5, where the hole could move only by its coupling either to orbitons or to magnons, but not to both of them simultaneously. Besides, the formation of the spin-orbital polarons was discussed in the lightly doped cobaltates but the orbital degrees of freedom there were integrated out; see M. Daghofer, P. Horsch, and G. Khaliullin, Phys. Rev. Lett.
**96**, 216404 (2006). - A. M. Oleś, P. Horsch, L. F. Feiner, and G. Khaliullin, Phys. Rev. Lett.
**96**, 147205 (2006). - S. Miyasaka, T. Okuda, and Y. Tokura, Phys. Rev. Lett.
**85**, 5388 (2000). - It is only the unusual coexistence of the $C\text{\u2212}\text{AF}$ phase and the metallic phase that was to some extent explained using the classical double exchange model adopted to the ${t}_{2g}$ orbital symmetries; see: K. Wohlfeld and A. M. Oleś, Phys. Status Solidi B
**243**, 142 (2006). - J. Fujioka, S. Miyasaka, and Y. Tokura, Phys. Rev. B
**72**, 024460 (2005);**77**, 144402 (2008). - P. Horsch, A. M. Oleś, L. F. Feiner, and G. Khaliullin, Phys. Rev. Lett.
**100**, 167205 (2008). - M. Imada, A. Fujimori, and Y. Tokura, Rev. Mod. Phys.
**70**, 1039 (1998). - J. van den Brink, G. Khaliullin, and D. Khomskii, in
*Colossal Magnetoresistive Manganites*, edited by T. Chatterji (Kluwer, Dordrecht, 2004). - S. Ishihara, Phys. Rev. Lett.
**94**, 156408 (2005). - A. M. Oleś, G. Khaliullin, P. Horsch, and L. F. Feiner, Phys. Rev. B
**72**, 214431 (2005). - E. Dagotto, Rev. Mod. Phys.
**66**, 763 (1994); M. Brunner, F. F. Assaad, and A. Muramatsu, Phys. Rev. B**62**, 15480 (2000); M. Bejas, A. Greco, and A. Foussats, ibid.**73**, 245104 (2006). - A. B. Harris, T. Yildirim, A. Aharony, O. Entin-Wohlman, and I. Y. Korenblit, Phys. Rev. B
**69**, 035107 (2004). - J. Fujioka, S. Miyasaka, and Y. Tokura, Phys. Rev. Lett.
**97**, 196401 (2006); M. De Raychaudhury, E. Pavarini, and O. K. Andersen, ibid.**99**, 126402 (2007). - T. Mizokawa and A. Fujimori, Phys. Rev. B
**54**, 5368 (1996). - G. Khaliullin, P. Horsch, and A. M. Oleś, Phys. Rev. B
**70**, 195103 (2004). - K. A. Chao, J. Spałek, and A. M. Oleś, J. Phys. C
**10**, L271 (1977); Phys. Rev. B**18**, 3453 (1978). - H. Eskes and A. M. Oleś, Phys. Rev. Lett.
**73**, 1279 (1994); H. Eskes, A. M. Oleś, M. B. J. Meinders, and W. Stephan, Phys. Rev. B**50**, 17980 (1994). - J. Bała, A. M. Oleś, and J. Zaanen, Phys. Rev. B
**52**, 4597 (1995). - A. M. Oleś, P. Horsch, and G. Khaliullin, Phys. Rev. B
**75**, 184434 (2007). - A. Weiße, J. Loos, and H. Fehske, Phys. Rev. B
**64**, 054406 (2001); J. Bała, G. A. Sawatzky, A. M. Oleś, and A. Macridin, Phys. Rev. Lett.**87**, 067204 (2001); A. Weiße and H. Fehske, New J. Phys.**6**, 158 (2004). - We have verified analytically that these terms would not lead to any $\mathbf{k}$ dependence in the spectra as in the self-energies describing such processes one can shift the summation over the momenta in such a way that the self-energies are momentum independent (See also Sec. III B for a similar calculation concerning the self-energies for two boson lines). Thus, these terms would not change the qualitative feature of the calculated spectra, i.e., that the $\mathbf{k}$ dependence of the quasiparticle states originates entirely from the three-site terms. Besides, we note that inclusion of such terms would require going beyond the LSW and LOW approximation, i.e., we would have to include the interactions between magnons and/or orbitons.
- Note, that only for didactic purpose in Secs. V A and V B we calculate the spectral functions for $\eta =0$. However, one should bear in mind that such a value is unphysical and inconsistent with the assumed infinite value of the Hund’s coupling in Eq. (2.2).
- S. Miyasaka, J. Fujioka, M. Iwama, Y. Okimoto, and Y. Tokura, Phys. Rev. B
**73**, 224436 (2006). - The obtained value of the incoherent bandwidth $W=4t$ in the classical limit well agrees with the retraceable path approximation (RPA) formula $W=4{t}_{\text{eff}}\sqrt{l}$, where the effective hopping ${t}_{\text{eff}}=t/\sqrt{2}$ due to the double exchange and $l=2$ is the number of possible forward going steps in the classical spin-orbital model, cf. discussion in the end of Sec. VC in Ref. 13. Note also that the narrowing of the bandwidth in the quantum case is due to the fact that the effective number of forward going steps, which the hole can make so that the spins become misaligned (which is the essence of RPA), is slightly reduced. This is because some of the spins are already overturned due to the quantum spin fluctuations in the full quantum model.
- J. Zaanen, A. M. Oleś, and P. Horsch, Phys. Rev. B
**46**, 5798 (1992). - W. F. Brinkman and T. M. Rice, Phys. Rev. B
**2**, 1324 (1970). - F. Krüger, S. Kumar, J. Zaanen, and J. van den Brink, Phys. Rev. B
**79**, 054504 (2009).