# Fidelity susceptibility in the two-dimensional spin-orbit models

###### Abstract

We study the quantum phase transitions in the two-dimensional spin-orbit models in terms of fidelity susceptibility and reduced fidelity susceptibility. An order-to-order phase transition is identified by fidelity susceptibility in the two-dimensional Heisenberg XXZ model with Dzyaloshinsky-Moriya interaction on a square lattice. The finite size scaling of fidelity susceptibility shows a power-law divergence at criticality, which indicates the quantum phase transition is of second order. Two distinct types of quantum phase transitions are witnessed by fidelity susceptibility in Kitaev-Heisenberg model on a hexagonal lattice. We exploit the symmetry of two-dimensional quantum compass model, and obtain a simple analytic expression of reduced fidelity susceptibility. Compared with the derivative of ground-state energy, the fidelity susceptibility is a bit more sensitive to phase transition. The violation of power-law behavior for the scaling of reduced fidelity susceptibility at criticality suggests that the quantum phase transition belongs to a first-order transition. We conclude that fidelity susceptibility and reduced fidelity susceptibility show great advantage to characterize diverse quantum phase transitions in spin-orbit models.

###### pacs:

03.67.-a,64.70.Tg,75.25.Dk,05.70.Jk## I Introduction

Spin-spin interactions have been intensively studied in quantum magnets and Mott insulators in the last decades Zaanen-PRL-55-418 (). The effect of the orbital degree of freedom has received much attention since the discovery of a variety of novel physical phenomena and a diversity of new phases in transition metal oxides (TMOs) Science288.462 (); NatMater6.927 (). In particular, under an octahedral environment, the orbital degeneracy of transition metal ions is partially lifted, and the remaining orbital degrees of freedom can be generally described by localized pseudospins. A typical effect induced by such a symmetry degradation of orbital degeneracy is the presence of bond-selective pseudospin interaction. It is because the spatial orientations of the orbits lead to anisotropic overlaps between neighboring ions. Consequently, the interactions among different bonds are intrinsically frustrated. To understand the orbital degree of freedom, the so-called quantum compass model (QCM) has been given rise to intensive research SovPhysJETP.37.725 (); PhysRevLett.93.207201 (); PhysRevB.71.024505 (); PhysRevB.72.024448 (); PhysRevB.71.195120 (); PhysRevB.75.144401 (); PhysRevB.75.134415 (); PhysRevB.78.064402 (); PhysRevLett.102.077203 (); PhysRevB.80.014405 (); PhysRevB.78.184406 (); PhysRevB.79.104429 (). In QCM, the pseudospin operators are coupled in such a way as to mimic the built-in competition between the orbital orders in different directions Science288.462 (); NatMater6.927 (). The frustration leads to macroscopic degeneracy in the classical ground state PhysRevLett.93.207201 (), and even highly degenerate quantum ground state PhysRevB.75.134415 (); PhysRevB.78.184406 (). Interestingly, the two-dimensional (2D) QCM has become a prototype to generate topologically protected qubits PhysRevB.71.024505 (); PhysRevLett.99.020503 ().

Besides, considering the spin-orbit coupling in TMOs is inevitable and intriguing. For example, octahedra tilt may give rise to effective Dzyaloshinsky-Moriya interaction (DMI) PhysRevLett.102.017205 (). Comparing those that only possess spin-spin coupling, spin-orbit models appear to be more intricate. The interplays between spin and orbital degree of freedom host a variety of different phases. For instance, the orbital exchange in a honeycomb lattice induces orbital ordering Wu-PRL-2006 (); PhysRevB.75.195118 () and topological order Kitaev ().

In this paper, we concentrate on portraying the quantum phase transitions (QPTs) in 2D spin-orbit interaction Hamiltonians. As we know, a QPT identifies any point of nonanalyticity in the ground-state (GS) energy of an infinite lattice system Sachdev (). Conventionally, local order parameters are needed to detect the nonanalyticity in the GS properties as the system varies across the quantum critical point (QCP). However, the knowledge of the local order parameter is not easy to retrieve from a general many-body system, especially for QPTs beyond the framework of the Landau-Ginzburg symmetry breaking paradigm Kitaev (); Feng-PRL-98-087204 (). Recently, quantum fidelity, also referred to as the GS fidelity, sparked great interest among the community to use it as a probe for the QCP HTQuan2006 (); Pzanardi2006 (). The fidelity defines the overlap between two neighboring ground states of a quantum Hamiltonian in the parameter space, i.e.,

(1) |

where is the GS wave function of a many-body Hamiltonian = + , is the external driving parameter, and is a tiny variation of the external parameter. Though borrowed from the quantum information theory, fidelity has been proved to be a useful and powerful tool to detect and characterize QPTs in condensed matter physics IntJModPhys.24.4371 (). In order to remove the artificial variation of external parameters, one of the authors and collaborators in Ref. PhysRevE.76.022101 () introduced the concept of fidelity susceptibility (FS), which is the leading-order term of fidelity with respect to the external driving parameter,

(2) |

FS elucidates the rate of change of fidelity under an infinitesimal variation of the driving parameter. There exists an intrinsic relation between the FS and the derivatives of GS energy as following:

(3) | |||||

(4) |

Here the eigenstates satisfy . The apparent similarity between Eq. (3) and Eq. (4) arouses similar critical behavior around the critical point.

In addition, the reduced fidelity and its susceptibility were also suggested in the studies of QPTs Paunkovic (); Jian-PRE-78-051126 (); Ho-PRA-78-062302 (); You-PLA (). The reduced fidelity concerns the similarity of a local region of the system with respect to the driving parameter. Despite the locality, the reduce fidelity fidelity (RFS) encodes the fingerprint of QPTs of the whole system, and is even more sensitive than its global counterparts Eriksson-PRA-2009 (); Zhi-PRA-2010 ().

For a second-order QPT, around the critical point, the correlation length diverges as , while the gap in the excitation spectrum vanishes as , where is the critical point in the thermodynamic limit, and are the critical exponents. At the critical point, i.e., , the only length scale is the system size . We can account for the divergence at the QPT by formulating a finite-size scaling (FSS) theory. Universal information could be decoded from the scaling behavior of FS Kwok-PRB-2008 (); Gu-PRB-2008 (); Yu-PRE-2009 (). FS increases as the system size grows, and the summation in Eq.(3) contributes to an extensive scaling of in the off-critical region. Therefore, FS per site appears to be a well-defined value, where is the number of sites and stands for system dimensionality. Instead, FS exhibits stronger dependence on system size across the critical point than in non-critical region, showing that a singularity emerges in the summation of Eq.(3). This implies an abrupt change in the ground state of the system at the QCP in the thermodynamic limit. Following standard arguments in scaling analysis Continentino (), one obtains that FS per site scales as Schwandt-PRL-2009 (); Grandi (); arXiv1104.4104 ()

(5) |

Due to the arbitrariness of relevance of the driving Hamiltonian under the renormalization group transformation, could be (i) superextensive if Schwandt-PRL-2009 (), (ii) extensive if Garnerone (); Jacobson-PRB-2009 (), (iii) subextensive if Venuti-PRL-99-095701 (). From Eq. (3) we can deduce that the ground state of the system should be gapless if is superextensive, but not vise versa. For finite size system, the position of a divergence peak defines a pseudocritical point as the precursor of a QPT, and it approaches the critical point as . For sufficiently relevant perturbations on large size system, i.e., , the leading term in expansion of pseudocritical point obeys such scaling behavior as Zanardi-PRA-2008 ()

(6) |

Thus, the behavior of on finite systems in the vicinity of a second-order QCP can be estimated as Schwandt-PRB-2010 (); arXiv1106.4078 (); arXiv1106.4031 ()

(7) |

where is an unknown regular scaling function, and is a constant independent of and .

However, many important QPTs fall in the category of first-order QPTs, in which the first derivative of GS energy exhibits discontinuity at critical point. Different from second-order QPT, there is no characteristic correlation length in first-order QPT, and in general the FSS will violate the scaling relations of second-order QPTs Binder-RPP-1987 (). For a typical first-order QPT, two competing ground states and are degenerate at critical point in the thermodynamic limit, and they become energetically favorable on one side of such as and respectively. The level crossings in the thermodynamic limit usually turn into avoided level crossings for finite system, and degeneracy at critical point in general is lifted, opening an energy gap . In the low-energy subspaces spanned by and , the diagonal matrix elements coincide =, but the off-diagonal matrix elements = induce an avoided level crossing with

(8) |

For a Hamiltonian with local interactions (e.g., nearest neighbors),

(9) |

where , , , () form a complete set of basis vectors in Hilbert space , and () are the corresponding coefficients. The dimension of is of ( is the degree of freedom of Hamiltonian constituent, e.g., for spin-1/2 Hamiltonian), and () are of order . The finite Hamming distance between and gives limited number (roughly speaking, ) of nonzero . Hence, the off-diagonal matrix elements scale exponentially , and consequently Schutzhold (); Jorg (). If any of the ground states are degenerate protected by symmetries, should be written in the subspace of all the degenerate states, and the dimension of the matrix will be larger than . Considering Eq. (3), the exponentially closing gap gives a dominant contribution, manifesting that the FS should carry the signature of an exponential divergence at critical point. Moreover, from another point of view, due to their macroscopic distinguishability, the overlap between states should be exponentially small with system size at criticality Schutzhold (). In general, FS per site scales exponentially, i.e.,

(10) |

where is a size independent constant, and is a polynomial function of , which is a correction to the exponential term.

In this work, we use FS and RFS to incarnate the critical phenomena in 2D spin-orbit models. We show that the FS and RFS manifest themselves as extreme points at the QCPs. In other words, FS and RFS are quite sensitive detectors of the QPTs in spin-orbit models. This article is organized as follows: In Sec. II, a general Hamiltonian is given. We investigate QPTs in the spin-orbit model on a square lattice in Sec. III. The phase diagram and the FSS are studied through FS. Next, in Sec. IV, we study QPTs in Kitaev-Heisenberg model on a 2D honeycomb lattice by FS and the second derivative of the the GS energy. Both approaches identify three distinct phases, and FS behaves more sensitively. After that, we take advantage of the RFS to locate the QCP of 2D compass model in Sec. V, and show that the exponentially divergent peaks of FS and RFS imply that the first order QPT. A brief summary is presented in Sec. VI.

## Ii General Hamiltonian

Motivated by a broad interest in Mott insulators, we construct the general structure of model Hamiltonian on a given bond connecting two nearest-neighbor (NN) sites and ,

(11) |

The first term in Hamiltonian (11) describes bond-selective interaction, where the orbital exchange constant is derived from multi-orbital Hubbard Hamiltonian consisting of the local on-site repulsion and the hopping integral . The second term corresponds to the isotropic Heisenberg coupling, and Dzyaloshinsky-Moriya anisotropy in the third term comes from the lattice distortions. The Hamiltonian (11) extrapolates from the Heisenberg model to QCM depending on the bonds geometry. In this paper, we devote our study to 2D lattices, either square or hexagonal lattice. Several real systems motivate the investigation of these kinds of lattices, such as SrTiO Schrimer-JPC-17-1321 (), SrIrO Cao-RPB-57-R11039 (); Moon-PRB-74-113104 (); Kim-PRL-101-076402 (). The flexibility of parameters induces a rich variety of phases and various fascinating physical phenomena.

## Iii Fidelity Susceptibility in the Spin-Orbit model

In strong spin-orbit coupling materials, an XXZ model with DMI was suggested,

(12) | |||||

where is the anisotropic parameter and is the strength of DMI. A general boundary condition can be written as , where , corresponds to open boundary condition (OBC) and to the periodic boundary condition (PBC). The Hamiltonian (12) was proposed to describe the layered compound SrIrO PhysRevLett.102.017205 (). The spin canting is induced by lattice distortion of corner-shared IrO octahedra. The DMI was introduced originally to explain the presence of weak ferromagnetism in antiferromagnetic (AFM) materials Dzyaloshinsky-JPCS-4-241 (); Moriya-PR-4-288 (), such as -FeO, MnCO and CoCO, since such antisymmetric interaction could produce small spin cantings. Recently, the influence of DMI has become very important in elucidating many interesting properties of different systems, e.g., ferroelectric polarization in multiferroic materials Sergienko-PRB-2006 (), exchange bias effects in perovskites Dong-PRL-2009 (), asymmetric spin-wave dispersion in double layer Fe Zakeri-PRL-2010 (); Michael-PRB-2010 (), and noncollinear magnetism in FePt alloy films Honolka-PRL-2009 (). For simplicity, we assume the DMI fluctuation is acting on the - plane, and the vector is imposed along the direction, i.e., .

The Hamiltonian (12) reduces to the anisotropic Heisenberg XXZ model when the rotations of IrO octahedra are absent, i.e., =0. When , one-dimensional (1D) spin-1/2 XXZ chain has long range order and gapped domain-wall excitations. On the other hand, in the XX limit, i.e., = 0, the system is equivalent to a chain of non-interacting fermionic model, which becomes gapless in the thermodynamic limit. The QPT takes place at the isotropic point = 1, which is a Berezinskii-Kosterlitz-Thouless (BKT) transition point. BKT phase transition belongs to an infinite-order phase transition, and the ground-state energy and all of its derivatives with respect to are continuous at the critical point. However, the FS succeeds in detecting the nonanalyticities of the ground state across BKT transition Yang-PRB-2007 (). For instance, the BKT transition from gapless Tomonaga-Luttinger liquid to gapped Ising phase in 1D XXZ model is detected by the divergence of the FS using the density-matrix-renormalization-group (DMRG) technique Wang-PRA-2010 (). In addition, the BKT transition from spin fluid to dimerized phase in the - model Chen-PRA-2008 () and the superfluid-insulator transition in the Bose-Hubbard model at integer filling Gu-PRB-2008 () are also able to be signaled by FS.

Different from the 1D case, 2D XXZ model exhibits a second-order QPT at the isotropic point =1 in the thermodynamic limit, where the first excited energy levels cross Tian-PRB-2003 (). For , the Ising term in the Hamiltonian dominates and the ground state is an AFM phase along the direction. For , the first two terms in the Hamiltonian dominate and the ground state is also an AFM phase, but in the - plane. It is well known that long-range orders are present in both phases Kennedy-PRL-1988 (). FS can serve as a sensitive detector of the critical point in the 2D XXZ model on a relatively small square lattice Yu-PRE-2009 ().

Actually, the DMI does not change the universality class of XXZ model. The DMI can be eliminated from the Hamiltonian (12) by a spin axes rotation Perk (); Aristov-PRB-2000 (); Alcaraz-JSP-1990 (), and after rotation an unitarily equivalent form is given by

(13) |

where , and with boundary condition . Note that the mapping becomes exactly equivalent only in the thermodynamic limit and for open boundary condition. In the thermodynamic limit, the boundary condition does not affect the critical behavior and consequently the will have the same critical properties as the . Hence, QCP at =1 in becomes a critical line in .

FS and the second derivative of GS energy of square lattice (see Fig. 1) as a function of anisotropy for and are demonstrated in Fig. 2. As is sketched, both approaches display broad peaks at an identical . However, the peak of FS is slightly narrower than that of second derivative of GS energy. The locations of peaks are specified as pseudocritical points. The heights of local maxima decrease as increases.

Next, to check whether these extreme points can be regarded as QCPs, we should perform FSS analysis, as is illustrated in Fig. 3. In Figs. 3(a) and 3(b), we plot FS per site as a function of for different lattice size =10, 16, 18, 20 and 26. As is clearly shown that becomes more pronounced for increasing . is extensive in the off-critical region, while superextensive at the pseudocritical point. From the scaling relation Eq. (5), a linear dependence for the maximum fidelity susceptibilities on is expected with effective length . This is confirmed by the results shown in Figs. 3(c) and 3(d), in which we plot the logarithm of maximum fidelity susceptibilities at pseudocritical points verse a function of . For large , scales linearly with , and the slope is interpreted as the inverse of critical exponent of the correlation length from Eq. (5). By applying linear regression to the raw data obtained from XXZ model on various square lattices, we derive for , while for . The values got from the measurements of and are consistent with each other up to two digits. The validation of Eq. (5) implies that the QPT in 2D XXZ with DMI manifests itself as a clear sign of a second-order QPT.

From Figs. 3(a) and 3(b), the positions of cusps seemingly converge toward the critical points. However, these system sizes seem insufficiently to estimate the critical exponents from Eqs. (6) and (7), and the corrections to scaling relations are not negligible. In this case, the shift of the location of the pseudocritical point from real critical point should be replaced with , where and are the coefficients Zanardi-PRA-2008 (). As is shown in insets of Figs. 3(c) and 3(d), we notice that the data points corresponding to the small system sizes clearly deviate from the linear fit obtained for the points for the two largest . As , the pseudocritical points approach the critical points . The linear fits from results of and (see Fig. 1) yield the estimates the are and , respectively. To get more precise critical exponents, more elaborate schemes are needed, such as large-scale quantum Monte Carlo simulations Schwandt-PRL-2009 (); Schwandt-PRB-2010 ().

To proceed, we calculate the FS on =20 square lattice for different sets of parameters and . The FS as a function of for and are depicted in Fig. 4, and the FS as a function of for =0, 1, 2, and 3 are plotted in Fig. 5. A boost of and suppress the FS. By retrieving the critical points from the locations of FS peaks, we are able to draw the corresponding phase diagrams. The pseudocritical lines in the plane and the plane are shown in the insets of Fig. 4 and Fig. 5. We can notice that pseudocritical line for (solid dotted line) are qualitatively similar to the critical line in the thermodynamic limit (dash line).

## Iv Fidelity susceptibility in the 2D Kitaev-Heisenberg model

Since Kitaev introduced a spin 1/2 quantum lattice model with Abelian and non-Abelian topological phases, finding a physical realization of Kitaev model has triggered a tremendous amount of interest PhysRevLett.105.027204 (). It is proposed that in iridium oxides AIrO, the strong spin-orbit coupling may lead to the desired anisotropy of the spin interaction. The Ir ions in iridium oxides AIrO can be effectively illustrated as spin half on a honeycomb lattice, where three distinct types of NN bonds are referred to as (, , ) bonds. To describe the competition between direct exchange and superexchange mechanisms, a so-called Kitaev-Heisenberg model was dedicated, i.e.,

(14) |

The Hamiltonian (14) has been parameterized, which encompasses a few well-known models. For , it reduces to the Heisenberg model on a hexagonal lattice, while it becomes exactly solvable ferromagnetic Kitaev model at . Another solvable point corresponds to , where it is unitarily equivalent to ferromagnetic Heisenberg model. The ground state of Hamiltonian (14) evolves from 2D Néel AFM state () to stripy antiferromagnetism (), and to spin liquid () as increases. We perform exact diagonalization (ED) on two geometries with =16 and =24 (see Fig. 6) to calculate the GS energy and fidelity susceptibility .

The second derivative of energy density and FS per site as a function of coupling strength are obtained by Lanczos calculation. As is illustrated in Fig. 7, two anomalies emerge when increases from to , indicating the system undergoes two phase transitions. The characterization of QPTs by FS is compatible with the second derivative of energy density. We find that they play similar roles in identifying the QPTs Chen-PRA-2008 (). However, the former demonstrates more pronounced peaks than the latter (notice the -axis is of a logarithmic scale). QPT from Néel-ordered to stripe-ordered phase takes place at , and this order-to-order phase transition seems insensitive to system size for these two configurations. On the other hand, the second order-to-disorder phase transition at appears to be sensitive to . The numerical calculations point out that the QPT around is first-order, while there is no consensus on the character of the QPT around HCJiang (); Fabien (). However, it is difficult to carry on ED on larger cluster for the hexagonal geometry and substantiate the scaling behavior of 2D Kitaev-Heisenberg model unless sophisticated techniques are used Schwandt-PRL-2009 (); Schwandt-PRB-2010 ().

## V Reduced Fidelity Susceptibility in The 2D Compass model

In recent years, the 2D AFM QCM has attracted considerable attention due to its interdisciplinary character. On one hand, it plays an important role in describing orbital interactions in TMOs. On the other hand, 2D QCM is equivalent to Xu-Moore model Xu2004 () and toric code model in a transverse field Vidal2009 (), which could possibly be used to generate protected qubits realized by Josephson-coupled superconducting arrays PhysRevB.71.024505 ().

QCM is defined on a square lattice with PBC by the Hamiltonian

(15) |

where are unit vectors along the and directions, and () is the coupling in the () direction. (=, , ) are the pseudospin operators of lattice site obeying . There exists a first-order phase transition at the self-dual point, i.e., , and it is widely believed that the phase transition belongs to first-order PhysRevB.72.024448 (); PhysRevB.75.144401 (); PhysRevLett.102.077203 (); Vidal2009 ().

In the case of being even, this model is equivalent to the ferromagnetic QCM by rotating the pseudospin operators at one sublattice by an angle about -axis. The Hamiltonian (15) enjoys such commutation relation with the parity operators and , which are defined as

(16) |

where indices and are the - and -component of lattice site . However, column parity operator does not commute but anticommutes with row parity operator . In such circumstance, the Hilbert space can then be decomposed into subspace , where is the eigenvalue of and is specified in each subspace. Since , can be either 1 or -1. The lowest energy of the Hamiltonian in the subspace is nondegenerate. In Ref. Brzezicki-PRB-2010 (), the authors unraveled the hidden dimer order therein. The symmetries in these systems induce large degeneracies in their energy spectra, and make their numerical simulation tricky. Recently, the ground state of is proven to reside in the most homogeneous subspaces and , in other words, they are two-fold degenerate You-JPA-2010 ().

At this stage, a single-site pseudospin-flipping interaction will change the parity of the chain, contrary to the flipping terms in . For example, or acting on one chain along direction will change the parity of the chain. For two-point correlation functions between sites and separated by sites, it is not difficult to yield

(17) |

Similarly, =0, , . The only correlation functions surviving are and , as are shown in Fig. 8 and inset of Fig. 11. As increases, decreases while increases accordingly. They become equal at = .

With these properties of the ground state of the 2D QCM, a simple form for the two-site density matrix is obtained,

(18) |

in which the prime means tracing over all the other pseudospin degrees of freedom except the two sites and . are Pauli matrices , and for = 1 to 3, and 2 by 2 unit matrix for =0. If the two pseudospins are linked by an -type bond, we have

(19) |

in which and are the 2 by 2 unit matrices. If translational invariance is preserved in the ground state , the reduced density matrix can be simplified as

(20) |

Similarly, if the two sites and are linked by -type bond, we have

(21) |

In addition, the NN two-point correlation functions can be calculated using the Feymann-Hellmann Theorem

(22) |

where means -type and -type bond respectively.

The reduced fidelity is defined as the overlap between and , i.e.

(23) |

The prime means a different reduced density matrix induced by the tiny changes of driving parameters in the Hamiltonian. By Eqs. (20) and (21), it is obvious that commutes with , and we can easily evaluate the reduced fidelity

(24) | |||||

The two-site RFS could also be obtained straightforwardly,

(25) |

in which (=,) is the driving parameter for the QPT. According to Eqs. (22) and (25), we observe that the numerator of is proportional to the square of the second derivative of GS energy and the denominator is finite in unpolarized state. The second power in the numerator indicates that the two-site RFS is more effective than the second derivative of the GS energy in measuring QPTs Xiong-PRB-2009 (). We compare the RFS of two NN sites with the second derivative of GS energy of square lattice (see Fig. 1) displayed in Fig. 9, and find that they indeed present cusp-shaped peaks, but the former is more sensitive to phase transition than the latter. The pronounced maxima of peaks arise at .

Fig. 10 unveils the results from ED on various 2D square lattices with PBC. As the coupling changes, both FS and RFS of two NN sites exhibit peaks around . With increasing system size, the peaks of (or ) become more pronounced and pseudocritical points seemingly converge toward the real critical point quickly. We plot the maximum fidelity and maximum reduced fidelity susceptibility against the system size in Fig. 10(d). The scaling of the peaks reveals approximately exponential divergences at criticality instead of power-law divergences, which indicates the occurrence of non-second-order phase transition.

As is argued in Eqs. (8) - (10), the approximately exponential divergence of FS and RFS hints an exponentially small gap at the critical point and thus a first-order phase transition. For two next-nearest-neighbor (NNN) sites, we calculate the RFS according to Eq. (25). We observe a similar behavior of the RFS with respect to , and find that the RFS of NNN pair is a bit larger than the RFS of NN bond. This is illustrated in Fig. 11. In a word, the two-site RFS can serve as a signature for the QPTs in the 2D compass model.

## Vi Summary and Discussion

Motivated by recent theoretical and experimental work on orbital degree of freedom in Mott insulators, we have studied the QPTs in various 2D spin-orbit models using FS and RFS. The numerical analysis is performed on 2D clusters with the Lanczos algorithm. The spin-orbit model hosts plentiful phases, including phase transitions within and beyond the framework of Landau-Ginzburg paradigm. For the 2D XXZ model with DMI, with increasing driving parameters, the system undergoes a second-order phase transition from AFM state along - plane to AFM state along direction. We compare FS with second derivative of GS energy, and find both of them exhibit similar peaks. The FSS demonstrates that FS per site should diverge in the thermodynamic limit at pseudocritical point, and the locations of extreme points approach QCPs accordingly. The power-law divergence of FS at criticality indicates the quantum phase transition is of second order and the critical exponent is obtained. Analogously, as the exchange coupling changes, the ground state of 2D Kitaev-Heisenberg model evolves from Néel AFM state to stripy AFM state, and to Kitaev spin liquid. The QCPs could be signaled by the peaks of both FS and second derivative of GS energy. The nonlocal symmetries in 2D AFM QCM guarantee that the RFS can be written in an analytical form as a function of the correlation functions. Peaks of FS and two-site RFS take place around =. The quasicritical point develops into critical point quickly with increasing the system size. Scaling of the peaks reveals an exponential divergence at criticality, which suggests that a first-order phase transition happens. In conclusion, the FS and RFS are effective tools in detecting diverse QPTs in 2D spin-orbit models, and their scaling behaviors may hint the orders of phase transitions.

## Vii Acknowledgments

We acknowledge useful discussions with Wing-Chi Yu and J. Sirker. Wen-Long You acknowledges the support of the Natural Science Foundation of Jiangsu Province under Grant No. 10KJB140010 and the National Natural Science Foundation of China under Grant No. 11004144. Yu-Li Dong acknowledges the support of the Specialized Research Fund for the Doctoral Program of Higher Education (Grant No. 20103201120002), Special Funds of the National Natural Science Foundation of China (Grant No. 11047168) and the National Natural Science Foundation of China (Grant No. 11074184).

## References

- (1) J. Zaanen, G. A. Sawatzky, and J.W. Allen, Phys. Rev. Lett. 55, 418 (1985).
- (2) Y. Tokura, and N. Nagaosa, Science. 288, 462 (2000).
- (3) S. W. Cheong, Nature Mater. 6, 927 (2007).
- (4) Anup Mishra, Michael Ma, Fu-Chun Zhang, Siegfried Guertler, Lei-Han Tang and Shaolong Wan, Phys. Rev. Lett. 93, 207201(2004).
- (5) K. I. Kugel and D. I. Khomskii, Sov. Phys. -JETP 37, 725 (1973).
- (6) Z. Nussinov and E. Fradkin, Phys. Rev. B 71, 195120 (2005).
- (7) J. Dorier, F. Becca, and F. Mila, Phys. Rev. B 72, 024448 (2005).
- (8) H.-D. Chen, C. Fang, J.-P. Hu, and H. Yao, Phys. Rev. B 75, 144401 (2007).
- (9) Román Orús, Andrew C. Doherty and Guifré Vidal, Phys. Rev. Lett. 102, 077203 (2009).
- (10) Sandro Wenzel and Wolfhard Janke, Phys. Rev. B 78, 064402 (2008).
- (11) Wojciech Brzezicki, Andrzej M. Oleś, Phys. Rev. B 80, 014405 (2009).
- (12) Ke-Wei Sun, Yu-Yu Zhang and Qing-Hu Chen, Phys. Rev. B 79, 104429 (2009).
- (13) W. Brzezicki, J. Dziarmaga, and A. M. Oleś, Phys. Rev. B 75, 134415 (2007).
- (14) Wen-Long You and Guang-Shan Tian, Phys. Rev. B 78, 184406 (2008).
- (15) B. Douçot, M. V. Feigel’man, L. B. Ioffe, and A. S. Ioselevich, Phys. Rev. B 71, 024505 (2005).
- (16) P. Milman, W. Maineult, S. Guibal, L. Guidoni, B. Douot, L. Ioffe, and T. Coudreau, Phys. Rev. Lett. 99, 020503 (2007).
- (17) G. Jackeli and G. Khaliullin, Phys. Rev. Lett. 102, 017205 (2009).
- (18) Wen-Long You, Guang-Shan Tian and Hai-Qing Lin, Phys. Rev. B 75, 195118 (2007).
- (19) Congjun Wu, Phys. Rev. Lett. 100, 200406 (2008).
- (20) A. Kitaev, Ann. Phys. (N.Y.) 321, 2 (2006).
- (21) S. Sachdev, Quantum Phase Transitions, (Cambridge University Press, Cambridge, UK, 2000)
- (22) Xiao-Yong Feng, Guang-Ming Zhang, and Tao Xiang, Phys. Rev. Lett. 98, 087204 (2007).
- (23) H. T. Quan, Z. Song, X. F. Liu, P. Zanardi, and C. P. Sun, Phys. Rev. Lett. 96, 140604 (2006).
- (24) P. Zanardi and N. Paunkovic, Phys. Rev. E 74, 031123 (2006).
- (25) Shi-Jian Gu, Int. J. Mod. Phys. B 24, 4371 (2010).
- (26) Wen-Long You, Ying-Wai Li and Shi-Jian Gu, Phys. Rev. E 76, 022101 (2007).
- (27) N. Paunkovic, P. D. Sacramento, P. Nogueira, V. R. Vieira, V. K. Dugaev, Phys. Rev. A 77 052302 (2008).
- (28) J. Ma, L. Xu, H. N. Xiong and X. Wang, Phys. Rev. E 78, 051126 (2008).
- (29) H. M. Kwok, C. S. Ho, S. J. Gu, Phys. Rev. A 78, 062302 (2008).
- (30) Wen-Long You and Wen-Long Lu, Phys. Lett. A 373, 1444 (2009) .
- (31) Zhi Wang, Tianxing Ma, Shi-Jian Gu, and Hai-Qing Lin,Phys. Rev. A 81, 062350 (2010).
- (32) Erik Eriksson and Henrik Johannesson, Phys. Rev. A 79, 060301(R)(2009).
- (33) H. M. Kwok, W. Q. Ning, S. J. Gu and H. Q. Lin, Phys. Rev. E 78, 032103 (2008).
- (34) Wing-Chi Yu, Ho-Man Kwok, Junpeng Cao, and Shi-Jian Gu, Phys. Rev. E 80, 021108 (2009).
- (35) S. J. Gu, H. M. Kwok, W. Q. Ning and H. Q. Lin, Phys. Rev. B 77, 245109 (2008).
- (36) M. A. Continentino, Quantum Scaling in Many-Body Systems (World Scientific Publishing, Singapore, 2001).
- (37) David Schwandt, Fabien Alet, and Sylvain Capponi, Phys. Rev. Lett. 103, 170501 (2009).
- (38) C. De Grandi, V. Gritsev, and A. Polkovnikov, Phys. Rev. B 81, 012303 (2010); Phys. Rev. B 81, 224301 (2010).
- (39) Marek M. Rams, Bogdan Damski, arXiv:1104.4104 (unpublished).
- (40) Silvano Garnerone, N. T. Jacobson, Stephan Haas, and Paolo Zanardi, Phys. Rev. Lett. 102, 057205 (2009).
- (41) N. T. Jacobson, Silvano Garnerone, Stephan Haas, and Paolo Zanardi, Phys. Rev. B 79, 184427 (2009).
- (42) L. Campos Venuti and P. Zanardi, Phys. Rev. Lett. 99, 095701 (2007).
- (43) Paolo Zanardi, Matteo G. A. Paris, and Lorenzo CamposVenuti, Phys. Rev. A 78, 042105 (2008).
- (44) Albuquerque, A. Fabricio and Alet, Fabien and Sire, Clément and Capponi, Sylvain, Phys. Rev. B 81, 064418 (2010).
- (45) C. De Grandi, A. Polkovnikov, A. W. Sandvik, arXiv:1106.4078 (unpublished).
- (46) Michael Kolodrubetz, David Pekker, Bryan K. Clark, Krishnendu Sengupta, arXiv:1106.4031 (unpublished).
- (47) Kurt Binder, Rep. Prog. Phys. 50, 783 (1987).
- (48) R. Schützhold, J. Low temp. Phys. 153, 228 (2008).
- (49) T. Jörg, F. Krzakala, J.Kurchan, A.C. Maggs and J. Pujos, EPL, 89, 40004 (2010).
- (50) O. F. Schirmer, A. Forster, H. Hesse, M. Wohleeke and S. Kapphan, J. Phys. C: Solid. State Phys. 17, 1321 (1984).
- (51) G. Cao, J. Bolivar, S. McCall, J. E. Crow, and R. P. Guertin, Phys. Rev. B 57, R11039 (1998).
- (52) S. J. Moon, M. W. Kim, K. W. Kim, Y. S. Lee, J.-Y. Kim, J.-H. Park, B. J. Kim, S.-J. Oh, S. Nakatsuji, Y. Maeno, I. Nagai, S. I. Ikeda, G. Cao, and T. W. Noh, Phys. Rev. B 74, 113104 (2006).
- (53) B. J. Kim, Hosub Jin, S. J. Moon, J.-Y. Kim, B.-G. Park, C. S. Leem, Jaejun Yu, T. W. Noh, C. Kim, S.-J. Oh, J.-H. Park, V. Durairaj, G. Cao, and E. Rotenberg, Phys. Rev. Lett. 101, 076402 (2008).
- (54) I. Dzyaloshinsky, J. Phys. Chem. Solids 4, 241 (1958).
- (55) T. Moriya, Phys. Rev 4, 288 (1960).
- (56) I. A. Sergienko and E. Dagotto, Phys. Rev. B 73, 094434 (2006).
- (57) Shuai Dong, Kunihiko Yamauchi, Seiji Yunoki, Rong Yu, Shuhua Liang, Adriana Moreo, J.-M. Liu, Silvia Picozzi, and Elbio Dagotto, Phys. Rev. Lett. 103, 127201 (2009).
- (58) Kh. Zakeri, Y. Zhang, J. Prokop, T.-H. Chuang, N. Sakr, W. X. Tang, and J. Kirschner, Phys. Rev. Lett. 104, 137203 (2010).
- (59) Thomas Michael and Steffen Trimper, Phys. Rev. B 82, 052401 (2010).
- (60) J. Honolka, T. Y. Lee, K. Kuhnke, A. Enders, R. Skomski, S. Bornemann, S. Mankovsky, J. Minár, J. Staunton, H. Ebert, M. Hessler, K. Fauth, G. Schütz, A. Buchsbaum, M. Schmid, P. Varga, and K. Kern, Phys. Rev. Lett. 102, 067207 (2009).
- (61) M. F. Yang, Phys. Rev. B 76, 180403(R) (2007).
- (62) Bo Wang, Mang Feng, and Ze-Qian Chen, Phys. Rev. A 81, 064301 (2010).
- (63) S. Chen, L.Wang, Y. Hao, and Y. Wang, Phys. Rev. A 77, 032111 (2008).
- (64) G. S. Tian and H. Q. Lin, Phys. Rev. B 67, 245105 (2003).
- (65) T. Kennedy, E. H. Lieb and B. S. Shastry, Phys. Rev. Lett. 61, 2582 (1988).
- (66) J. H. H. Perk and H. W. Capel, Phys. Lett. A 58 115 (1976).
- (67) D. N. Aristov and S. V. Maleyev, Phys. Rev. B 62, 751(R) (2000).
- (68) F. C. Alcaraz and W. F. Wreszinski, J. Stat. Phys. 58, 45 (1990).
- (69) Jiří Chaloupka, George Jackeli, and Giniyat Khaliullin, Phys. Rev. Lett. 105, 027204 (2010).
- (70) Hong-Chen Jiang, Zheng-Cheng Gu, Xiao-Liang Qi, and Simon Trebst, Phys. Rev B, 83, 245104 (2011).
- (71) Fabien Trousselet, Giniyat Khaliullin and Peter Horsch, Phys. Rev B, 84, 054409 (2011).
- (72) C. Xu and J. E. Moore, Phys. Rev. Lett. 93, 047003 (2004).
- (73) Julien Vidal, Ronny Thomale, Kai Philip Schmidt, and Sébastien Dusuel, Phys. Rev. B 80, 081104(R) (2009).
- (74) Wojciech Brzezicki and Andrzej M. Oleś, Phys. Rev. B 82, 060401(R) 2010.
- (75) W.-L. You, G.-S. Tian, and H.-Q. Lin, J. Phys. A 43, 275001 (2010).
- (76) Heng-Na Xiong,Jian Ma, Zhe Sun, and Xiaoguang Wang, Phys. Rev. B 79, 174425 (2009).