Statistical and Dynamical Properties in the Energy
Spectra of XXZ Spin Chains
Kazue Kudo
Integrated Sciences, Graduate School of Humanities and Sciences,
Ochanomizu University
January, 2005

Abstract
Spectral properties of XXZ spin chains are studied intensively. The integrable XXZ spin
chain is one of the most important quantum integrable systems, and its Hamiltonian isgiven by
H=J

L


j=1

y
x
z
(Sjx Sj+1+ Sjy Sj+1
+ ∆Sjz Sj+1).

(1)

Here, L is the number of sites; S jα = 12 σjα and σjα (α = x, y, z) is the Pauli matrix on
the jth site. Equation (1) corresponds to the Heisenberg spin chain at ∆ = 1, to the
XY model at ∆ = 0, and to the Ising spin chain when ∆ → ∞. In this thesis, we study
on XXZ spin chains with next-nearest-neighbor (NNN) coupling and/or with a magnetic
ﬁeld as well as Eq. (1).
First, we discuss level crossings in the energy spectra of the integrable XXZ spin
chain. Some of the level crossings are counterexamples of a theorem by von Neumann andWigner: two energy levels with the same symmetry rarely cross each other in the ﬂow
of the spectrum. Multiple degeneracies appear at the particular crossing points where
the symmetry of the system is enhanced in the spectral ﬂow. The eigenvalues and theeigenvectors are obtained in an analytical form for 3-site and 4-site chains. The spectral
ﬂow of a 6-site chain is shown by numerical calculation.
Not only level crossings but also “level bifurcations” are found in the ﬂow of eigenvalues for the inhomogeneous transfer matrix of the XXZ spin chain. The “level bifurcations”
are branches such as bifurcation or multiple branches in the ﬂow of eigenvalues and notobserved for the homogeneous transfer matrix. The observations are obtained by numerical calculation. We also conﬁrm the appearance of some of the branches through the
Bethe ansatz method. The behavior of some branching points are analyzed in the ﬂow of
eigenvalues for a particular inhomogeneous case.
Level statistics is discussed for the XXZ spin chains with NNN coupling or with arandom magnetic ﬁeld. Evaluating the level-spacing distribution, the spectral rigidity andthe number variance, we conﬁrm the following conjecture: the level-spacing distribution
should be described by the Poisson distribution for an integrable model and by the Wigner
distribution for a nonintegrable model. The XXZ spin chains with the NNN coupling (J 2 )or with a random magnetic ﬁeld (h) are nonintegrable, so that mostly the level statistics
shows Wigner behavior. However, the level statistics even for a nonintegrable model
sometimes shows Poisson-like behavior because of the ﬁnite-size eﬀects. For the NNN
coupled XXZ spin chain, the behavior of level statistics comes close to the Poissonianaround J2 = 0 and/or ∆ = 1. As h/J becomes large, Poisson-like behavior appears forthe level statistics of the XXZ spin chain under random magnetic ﬁeld. When h/J  0,
however, Poisson-like behavior appears for small L.
In addition to the Poisson-like behavior, unexpected non-Wigner behavior of levelstatistics can appear due to mixed symmetry. When we discuss level statistics, we haveto consider the whole symmetry of a system and perform desymmetrization. The NNN
coupled XXZ chains are invariant under spin rotation around the z axis, translation,z= 0. If
reﬂection. However, we have to consider spin reversal symmetry only when S totzspin reversal symmetry is not considered for S tot = 0, the behavior of level statistics is
described by almost the numerical average of Poisson and Wigner behavior.
To observe the dynamical properties of the energy spectra, we investigate energy
diﬀusion for NNN coupled XXZ spin chains in the presence of a periodically oscillating
magnetic ﬁeld. Diﬀusion coeﬃcients (D) are found to obey the power law with respect
to both the ﬁeld strength (B0 ) and driving frequency (ω): D ∝ (B0 ω)β . We observethat β = 2 for “linear response regime” and β = 1 for “non-perturbative regime”. In thelinear response regime, the diﬀusion coeﬃcient can be calculated using the Kubo formula.
On the other hand, B0 ω is so large in the non-perturbative regime that the perturbation
theory fails. The diﬀusion coeﬃcient is smaller than that predicted by the Kubo formulaand the perturbation theory fails in the non-perturbative regime. The ranges of the linear
response and the non-perturbative regimes depend on J2 . When J2 is small, we observe
large-amplitude oscillations, which is attributed to a non-diﬀusive and ballistic nature of
the underlying energy diﬀusion.
In this way, we investigate the energy spectra of ﬁnite quantum XXZ spin chainsin several diﬀerent aspects. We hope that the results of the thesis will contribute to thedevelopment of the physics of low-dimensional systems and mesoscopic systems.

2

Contents
1

Introduction

3

2

Level Crossing
7
2.1 The eigenvalues and the eigenvectors of 3-site and 4-site XXZ spin chains7
2.2 The spectral ﬂow of the XXZ spin chain on 6 sites. . . . . . . . . . . . . . 13

3

Level Bifurcation
14
3.1 Transfer matrix . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
3.2 Eigenvalues of the transfer matrix obtained by numerical diagonalization . 15
3.3 Eigenvalues of the transfer matrix obtained through the Bethe ansatz method 19
3.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20

4Level Statistics
23
4.1 Preliminaries . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
4.2 Models and numerical procedure . . . . . . . . . . . . . . . . . . . . . . . 25
4.2.1 Models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
4.2.2 Desymmetrization . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
4.2.3 Unfolding . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
4.2.4 Level-spacing distribution, spectral rigidity and number variance . . 28
4.3 Integrable XXZ spin chain . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
4.4 Next-nearest-neighbor coupled XXZ spin chains . . . . . . . . . . . . . . . 31z
= 0 . . . . . . . . . . . . . . . . . . . . . . 31
4.4.1 Wigner behavior for Stot
4.4.2 How does non-Wigner behavior appear? . . . . . . . . . . . . . . . 33
4.4.3 Finite-size eﬀects on the level spacing distribution . . . . . . . . . . 37
4.4.4 Homogeneity of the characteristic behavior of level statistics throughout the spectrum . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
4.5 XXZ spin chains under random magnetic ﬁeld . . . . . . . . . . . . . . . . 40
4.5.1 Dependence on the anisotropy parameter ∆ . . . . . . . . . . . . . 40
4.5.2 Dependence on the magnetic ﬁeld . . . . . . . . . . . . . . . . . . . 42
4.5.3 Dependence on the system size . . . . . . . . . . . . . . . . . . . . 44

1

4.5.4
4.6

Homogeneity of the characteristic behavior of level statistics throughout the spectrum . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46

5Energy Diﬀusion
47
5.1 Preliminaries . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
5.2 What is energy diﬀusion like? . . . . . . . . . . . . . . . . . . . . . . . . . 49
5.3 Numerical procedure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
5.4 Time dependence of energy variance . . . . . . . . . . . . . . . . . . . . . 51
5.5 Diﬀusion coeﬃcients: dependence on ﬁeld strength and frequency . . . . . 52
5.6 Oscillation of energy diﬀusion in weakly-frustrated cases . . . . . . . . . . 55
5.7 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60

6Conclusions

61

A

Hamiltonian Matrix

64

B

Jordan-Wigner and Fourier Transformations

66

C

Spin reversal symmetry on momentum-based fermions

69

D Spin reversal operation on the vacuum state

71

E

73

Energy diﬀusion for a random matrix model2

Chapter 1
Introduction
A repulsion of two close levels is often discussed for various quantum systems through theperturbation theory [1, 2, 3]. Now we consider a Hamiltonian which is given by
H = H0 + V,

(1.1)

where H0 is the non-perturbative part and V is the perturbative part. When H0 has two
close eigenvalues ε 1 and ε2 , where the corresponding eigenstates are φ1 and φ2 , we canassume that the eigenstates of H should be given byψ = c 1 φ1 + c2 φ2 .

(1.2)

Then, the time-independent Schrödinger equation isHψ = Eψ,
where



(1.3)


H11 H12 .
H =
H21 H22(1.4)

The solutions of the equation is written as
1
1
(H11 − H22 )2 + 4|H12 |2,(H11 + H22 ) +2
2
1
1
(H11 + H22 ) −=(H11 − H22 )2 + 4|H12 |2 .
2
2

E1 =

(1.5)

E2

(1.6)

The energy levels are drawn in Fig. 1.1 as the functions of x = H11 − H22 . If the twolevels have diﬀerent symmetries, the level repulsion does not necessarily occur, since we
can choose such a basis that the non-diagonal elements H12 and H21 equal 0. If thetwo levels have the same set of symmetry quantum numbers, however, the non-diagonal
elements are generally nonzero and the level repulsion occurs.
According to a theorem by von Neumann and Wigner, two energy levels with thesame symmetry rarely cross each other in the ﬂow of the spectrum [4]. Here we assume3

Energy

H22

E1
2|H12|

0H11

x=H11-H22

E2

Figure 1.1: Level repulsion. (See text.)

E

x2

x1
Figure 1.2: Energy surfaces which have the same symmetry. x 1 and x2 are real parameters.

that the Hamiltonian matrix depends on two real parameters. Then the energy surfaces
which have the same symmetry near the crossing point are given by two circular corns
as shown in Fig. 1.2. We suppose the spectral ﬂow is given by changing a parameter x2,
where x1 is ﬁxed. Generally, a level crossing does not happen (Solid lines in Fig. 1.2).
When x1 is chosen to be through the top of the corns, however, the level crossing occurs(Dashed lines in Fig. 1.2). In other words, levels of the same symmetry can cross in the
spectral ﬂow, although the possibility is very small. When the level crossing occurs, thelevels are degenerate at the crossing. We call such a degeneracy an accidental degeneracy.
Some examples of accidental degeneracies are discussed, which give counterexamples tothe noncrossing rule [5, 6]. In particular, it is shown by Heilmann and Lieb that the onedimensional Hubbard model with the Hubbard coupling U provides a quantum mechanicalcounterexample to the noncrossing rule [5].
The above arguments indicate that spectral properties are relevant to symmetry of a4system. In this thesis, spectral properties are studied for quantum XXZ spin chains, which
can have extraordinary symmetry, i.e. sl 2 loop algebra symmetry. The XXZ Hamiltonian
on L sites is given by
H=J

L

j=1

y
x
z
(Sjx Sj+1+ Sjy Sj+1
+ ∆Sjz Sj+1
).

(1.7)

Here, S jα = 12 σjα, and σjα’s (α = x, y, z) are the Pauli matrices acting on the jth site
for j = 1, · · · , L. The XXZ Hamiltonian commutes with the sl 2 loop algebra when q
is a root of unity, where q is deﬁned by ∆ = (q + q −1)/2 [7]. The loop algebra is an
inﬁnite-dimensional Lie algebra and the dimensions of the largest degenerate eigenspaces
are given by N 2L/N , where q 2N = 1 for an integer N [8, 9]. Thus the spectral degeneracy
can increase exponentially with respect to L.
There is another motivation to study about XXZ spin chains. Quantum XXZ spin
chains are related to various important systems of quantum spin chains. An XXZ spin
chain corresponds to the Heisenberg spin chain at ∆ = 1, to the XY model at ∆ = 0, andto the Ising spin chain when ∆ → ∞. Furthermore, XXZ spin chains with next-nearestneighbor couplings are realized in quasi one-dimensional antiferromagnets: for example,
CsCoBr3 , CsCoCl3 [10], Cu(ampy)Br2 [11] and (N2H5)CuCl3 [12].
This thesis is based on 6 articles, i.e. Refs. [13, 14, 15, 16, 17, 18]. In Ref. [13],we have observed several level crossings associated with the sl 2 loop algebra symmetry
in the spectral ﬂow of the XXZ spin chain. Furthermore, novel spectral phenomena havebeen found in Ref. [14]. We have observed that branches appear in the spectral ﬂow
for the inhomogeneous transfer matrix of the XXZ spin chain. The spectral branches
are important, since an inhomogeneous transfer matrix plays a central role in the exact
calculation of thermodynamic quantities for the spin chain.
In Ref. [15], we have discussed the level-spacing distributions of next-nearest-neighbor(NNN) coupled XXZ spin chains to investigate the statistical behavior of their whole energy spectra. It is well known that a level-spacing distribution is described by the Poisson
distribution for an integrable model and by the Wigner distribution for a nonintegrable
model. Though the NNN coupled XXZ spin chains are nonintegrable models, robust nonWigner behavior has been seen. To ﬁnd possible distribution of another nonintegrable
XXZ model, i.e. an XXZ spin chain under random magnetic ﬁeld, in Ref. [16]. We haveshown how the level-spacing distributions depend on the lattice size, the anisotropy parameter and the magnetic ﬁeld. In Ref. [17], we have solved the unexpected non-Wigner
behavior of the level-spacing distributions of the NNN coupled XXZ spin chains. We havedemonstrated that mixed symmetry and ﬁnite eﬀects can lead to the unexpected behaviorof level statistics. We have also shown that the characteristic behavior of level statistics
does not depend on the spectral range for the XXZ spin chains.
We have investigated dynamics for the XXZ spin chains in Ref. [18]. The nature5

of energy diﬀusion has been analyzed for NNN coupled XXZ spin chains with an applied
periodically oscillating magnetic ﬁeld. We have evaluated the time dependence of energy
variance and shown how the diﬀusion coeﬃcients depend on the coupling constants, theanisotropy parameter, the magnetic ﬁeld and the driving frequency.
The organization of this thesis is the following. We show some counterexamples
to the noncrossing rule in Chapter 2. The spectral ﬂows are calculated for the XXZ
spin chains on 3, 4, and 6 sites. Diagonalizing the Hamiltonian, we observe several levelcrossings corresponding to the sl2 loop algebra symmetry. In Chapter 3, novel spectral
phenomena are shown in the spectral ﬂow of the XXZ transfer matrix. One spectralline may split into two spectral lines in the spectral ﬂow of the inhomogeneous XXZ
transfer matrix. Such branches appear in the spectral ﬂow for the inhomogeneous case,while no branches appear for the homogeneous case. In Chapter 4, level statistics is
discussed for XXZ spin chains with next-nearest-neighbor couplings or with a randommagnetic ﬁeld. Level statistics is described by Poisson behavior for integrable models
and Wigner behavior for nonintegrable ones. We investigate the correspondence between
non-integrability and Wigner behavior. We also refer that level statistics is inﬂuenced byﬁnite-size eﬀects and the symmetry of the system. Dynamical properties of the energy
spectra are discussed in Chapter 5. As dynamics of energy levels, energy diﬀusion isinvestigated. Diﬀusion coeﬃcients are found to obey the power law with respect to boththe ﬁeld strength and driving frequency with its power varying depending on the linear
response and non-perturbative regimes. Finally, we give conclusions in Chapter 6.

6

Chapter 2
Level Crossing
Two energy levels with the same symmetry rarely cross each other. The integrable XXZ
spin chain, however, has many level crossings at some particular values of the anisotropy
parameter associated with roots of unity. Diagonalizing the Hamiltonian of the XXZ spin
chain, we indeed observe some of such level crossings [13].
In this chapter, the eigenvalues and the eigenvectors are resolved analytically for
3-site and 4-site XXZ spin chains under the twisted boundary conditions. It will be usefulin the future to show analytical forms of the eigenvalues and the eigenvectors of XXZ
spin chains, even though the size of the chains is small. Furthermore, the spectral ﬂow isshown for the 6-site XXZ spin chain. Then several level crossings related to the sl 2 loop
algebra symmetry are observed.

2.1

The eigenvalues and the eigenvectors of 3-site
and 4-site XXZ spin chains

First of all, let us rewrite the Hamiltonian (1.7) by using σ j+ = (σjx + iσjy )/2 and σj− =(σjx − iσjy )/2:LJ+z
(2σ + σ − + 2σj− σj+1
+ ∆σjz σj+1
).
(2.1)
H=
4 j=1 j j+1±zz= e±iφ S1± , SN
Here, twisted boundary conditions are imposed: S L+1+1 = S1 . The conditions correspond to the state where a magnetic ﬂux φ goes through the spin chain. The
schematic picture of this model is shown in Fig. 2.1. In Appendix A, we show how towrite the Hamiltonian matrix of Eq. (2.1).
The eigenvalues and eigenvectors are resolved analytically. They are listed in Tables 2.1 and 2.2 for L = 3 and L = 4, respectively. The eigenvectors are not normalized.
Figures 2.2 and 2.3 shows the energy spectra for L = 3 and L = 4, respectively.
All the eigenvalues for L = 3 are double degenerated. The degeneracies correspond

7

5

φ

4

1

2

N

3

Figure 2.1: Schematic picture of the XXZ spin chain described by Eq. (2.1).

Table 2.1: The eigenvectors and the eigenvalues for the 3-site XXZ chain. The numbersin the ﬁrst column “No.” correspond to the numbers in Fig. 2.2.
No.
1
2

Eigenvectors| ↑↑↑| ↓↓↓

34

e 3 | ↑↑↓ + | ↑↓↑ + e− 3 | ↓↑↑φiφie 3 | ↑↓↓ + | ↓↑↓ + e− 3 | ↓↓↑

5
6
7
8

φi

φ−π

Eigenvalues(J/4) · 3∆
(double degenerated)
φiφ−π(J/4){4 cos[φ/3] − ∆}(double degenerated)e 3 i | ↑↑↓ − | ↑↓↑ + e− 3 i | ↓↑↑ (J/4){−4 cos[(φ − π)/3] − ∆}φ−πφ−πe 3 i | ↑↓↓ − | ↓↑↓ + e− 3 i | ↓↓↑(double degenerated)
φ+πφ+πe 3 i | ↑↑↓ − | ↑↓↑ + e− 3 i | ↓↑↑ (J/4){−4 cos[(φ + π)/3] − ∆}φ+πφ+πe 3 i | ↑↓↓ − | ↓↑↓ + e− 3 i | ↓↓↑(double degenerated)

8

1
0.8
3,4
0.6
1,2
0.4
0.2
E 0
-0.2
-0.4 7,8
-0.6
5,6
-0.8
-1
-1.2
0
1
0.8
3,4
0.6
0.4 1,2
0.2
E 0
-0.2
7,8
-0.4
-0.6
5,6
-0.8
-1
-1.2
0

1.5

(a)

(b)

1 3,4
0.5
E

0

5-8

-0.5

π/2

φ

π

3π/2

-1
-1

2π

1,2
-0.5

∆

0

0.5

1

1.5
3,4
1

(c)

(d)

0.5 7,8
E

0

-0.5
-1
π/2

φ

π

3π/2

2π

-1.5
-1

1,2
5,6

-0.5

∆

0

0.5

1

Figure 2.2: Energy spectra of XXZ spin chains on 3 sites. The numbers in the ﬁguresz|: Red dashed lines,
correspond to “No.” of Table 2.1. Color of the lines indicates |Stotzz| = 3/2; Blue solid lines, |S tot| = 1/2. The values of the ﬁxed parameters are the|Stotfollowing: (a) ∆ = cos(π/3), (b) φ = 0, (c) ∆ = 0.3, (d) φ = 0.6π.

9

1.5
1

2

15

1,2

0.5
5,6
E

13

0.5
E

3,4
14

14
1,2

7,8

-1

7,8
-1.5 16

16

-1.5
π/2

φ

π

3π/2

-2
-1

2π

1.5

-0.5

0.5

1

13
1,2

1

(c)

5,6

9,10

0.5
3,4

(d)

13

9,10

0

11,12

E

-0.5 14

5,6

0

-1

7,8

1,2

-1.5

11,12

14
3,4

-0.5

-1

-2
0

0

15

1.5

0.5
E

∆

2
15

1

13

3-6,11,12

0

-0.5

-1

-2
0

(b)

1

(a)

11,12

0

-0.5

15
9,10

1.5
9,10

7,8
16

-1.5
16
π/2

φ

π

3π/2

2π

-2
-1

-0.5

∆

0

0.5

1

Figure 2.3: Energy spectra of XXZ spin chains on 4 sites. The numbers in the ﬁguresz|: Green chained lines,
correspond to “No.” of Table 2.2. Color of the lines indicates |Stotzzz| = 2; Red dashed lines, |S tot| = 1; Blue solid lines, S tot= 0. The values of the ﬁxed
|Stotparameters are the following: (a) ∆ = cos(π/4), (b) φ = 0, (c) ∆ = 0.3, (d) φ = 0.6π.

to the Kramers doublets. Examples of the energy spectral ﬂows for L = 3 are shown in
Fig. 2.2. The left and right ﬁgures are drawn for the ﬁxed values of ∆ and φ, respectively.
In Fig. 2.2(a), where ∆ = cos(π/3), a sixfold degenerated point exists at φ = π. In
Fig. 2.2(b), where φ = 0, a sixfold degenerated point exists at ∆ = −0.5 = − cos(π/3).
There are not sixfold degenerated points for general cases, whose examples are shown in
Figs. 2.2(c) and 2.2(d). From Fig. 2.2 and the Table 2.1, we ﬁnd that sixfold degenerated
points appear for (∆, φ) = (− cos π3 , 0) and (cos π3 , π). We note that multi-fold degeneracies
appear when the system has many symmetries. Figures 2.2(a) and 2.2(b) are for the cases
where the system has many symmetries, while Figs. 2.2(c) and 2.2(d) are for the cases
where the system has fewer symmetries.
zz| = 2 and |Stot| = 1 are double degenerated
When L = 4, all the eigenvalues for |S totz= 0. It is a nature of even-site spin chains under zero magnetic ﬁeld.
but not for Stot
Examples of the energy spectral ﬂows for L = 4 are shown in Fig. 2.3. The left and10Table 2.2: The eigenvectors and the eigenvalues for the 4-site XXZ chains. The numbersin the ﬁrst column “No.” correspond to the numbers in Fig. 2.3.
No.
1
2

Eigenvectors| ↑↑↑↑| ↓↓↓↓3

eφi | ↑↑↑↓ + e( 4 φ− 2 )i | ↑↑↓↑ − e 2 | ↑↓↑↑ + e( 4 + 2 )i | ↓↑↑↑
φ π
3
π
φi
eφi | ↑↓↓↓ + e( 4 φ− 2 )i | ↓↑↓↓ − e 2 | ↓↓↑↓ + e( 4 + 2 )i | ↓↓↓↑

4
5
6
7
8

EigenvaluesJ∆(double degenerated)

π

φi

φ

3
φ+ π2
4

φi
2

φ π
−2
4

3

π

)i | ↑↑↓↑ − e | ↑↓↑↑ + e(
)i | ↓↑↑↑
eφi | ↑↑↑↓ + e(
φ π
3
π
φi
eφi | ↑↓↓↓ + e( 4 φ+ 2 )i | ↓↑↓↓ − e 2 | ↓↓↑↓ + e( 4 − 2 )i | ↓↓↓↑(J/2)(−2 sin(φ/4))(double degenerated)(J/2)(2 sin(φ/4))

e | ↑↑↑↓ − e | ↑↑↓↑ + e | ↑↓↑↑ − e | ↓↑↑↑
3φi
φi
φi
eφi | ↑↓↓↓ − e 4 | ↓↑↓↓ + e 2 | ↓↓↑↓ − e 4 | ↓↓↓↑

(double degenerated)(J/2)(−2 cos(φ/4))(double degenerated)

9
10
11
12

eφi | ↑↑↑↓ + e 4 | ↑↑↓↑ + e 2 | ↑↓↑↑ + e 4 | ↓↑↑↑
3φi
φi
φi
eφi | ↑↓↓↓ + e 4 | ↓↑↓↓ + e 2 | ↓↓↑↓ + e 4 | ↓↓↓↑| ↑↓↓↑ − | ↓↑↑↓
φiφie 2 | ↑↑↓↓ − e− 2 | ↓↓↑↑




φi
φi
φi
eφi − e 2 | ↑↑↓↓ − λ13 e 2 | ↑↓↑↓ + 1 − e 2 | ↑↓↓↑(J/2)(2 cos(φ/4))(double degenerated)0(double degenerated)(J/2) {−∆ + [ 413+ 1 − e 2 | ↓↑↑↓ + λ13 | ↓↑↓↑ + 1 − e− 2 | ↓↓↑↑

φi



φi

3φi
4

φi
2

φi
4

3φi

φi

φi





(Here, λ13 = −∆ +




φi

φi





φi

4 − 4 cos φ2 + ∆2.)


φi



eφi − e 2 | ↑↑↓↓ − λ14 e 2 | ↑↓↑↓ + 1 − e 2 | ↑↓↓↑

14



+ 1−e


15

16





| ↓↑↑↓ + λ14 | ↓↑↓↑ + 1 − e(Here, λ14 = −∆ −eφi + e


φi
2

+ 1+e


φi
2



φi
2

eφi + e


+ 1+e



φi
2



| ↓↓↑↑



φi





| ↓↑↑↓ + λ15 | ↓↑↓↑ + 1 + e


− φi
2



| ↓↓↑↑| ↓↑↑↓ + λ16 | ↓↑↓↑ + 1 + e(Here, λ16 = −∆ −1/2

−4 cos(φ/2) + ∆2 ](J/2) {−∆ + [ 4+4 cos(φ/2) + ∆2]

1/2

4 + 4 cos(φ/2) + ∆2 .)

φi
2

φi



| ↑↑↓↓ + λ16 e | ↑↓↑↓ + 1 + e 2 | ↑↓↓↑(J/2) {−∆ − [ 44 − 4 cos(φ/2) + ∆2 .)φi
2(Here, λ15 = −∆ +
φi
2



− φi
2

| ↑↑↓↓ + λ15 e | ↑↓↑↓ + 1 + e 2 | ↑↓↓↑

1/2−4 cos(φ/2) + ∆2 ]



− φi
2



| ↓↓↑↑4 + 4 cos(φ/2) + ∆2.)11(J/2) {−∆ − [ 4+4 cos(φ/2) + ∆2]

1/2

right ﬁgures are drawn for the ﬁxed values of ∆ and φ, respectively. In Fig. 2.3(a), where∆ = cos(π/4), an eightfold degenerated point exists at φ = π. For general values of ∆,the eightfold degenerated point vanishes and a pair of ﬁvefold degenerated points appears[see Fig. 2.3(c)]. The coordinates of the degenerated points depend on the value of ∆. In
Fig. 2.3(b), where φ = 0, a tenfold and two ﬁvefold degenerated points exist at ∆ = 0 andat ∆ = ±1, respectively. For general values of φ, the ten fold degenerated point vanishes
and an additional pair of ﬁvefold degenerated points appears [see Fig. 2.3(d)].
In this way, we ﬁnd many degeneracies at several points where the parameters, ∆ andφ, correspond to some particular values. Some of the particular values are independent of
the lattice size L (for example, φ = 0 or ∆ = 0). Other are dependent on L (for example,∆ = cos(π/3) or ∆ = cos(π/4)).
The analytical forms of the eigenvalues are useful to discuss the behavior of crossing
points in the energy spectra. In Fig. 2.3, for example, Line 1,2, Line 5,6 and Line 13 seemto be always crossing together at one point. We prove that the expectation is right from
Table 2.2. At the crossing point of Line 1,2 and Line 5,6,φ
∆ = sin .
4(2.2)

Here, we assume J = 0. On the other hand, at the crossing point of Line 1,2 and Line 13,


1 

2
1
φ
∆=
−∆ + 4 − 4 cos + ∆2
.

2
2

Namely,

(2.3)

1

2
φ
3∆ = 4 − 4 cos + ∆2 .
2

(2.4)

Here we notice that ∆ ≥ 0. If each hand side of Eq. (2.4) is squared, then we have


∆

2

1
φ
=
1 − cos22
φ= sin2 .
4



(2.5)
(2.6)

Since ∆ ≥ 0, Eq. (2.6) is reduced toφ∆ = sin .
4(2.7)

We now ﬁnd that Eq. (2.2) coincides with Eq. (2.7). In other words, the crossing point ofLine 1,2 and Line 5,6 coincides with that of Line 1,2 and Line 13. Therefore those lines
are always crossing together at one point.

12

3
2
E 1
0
-1
-2
-3
-1

-0.5

∆

0

0.5

1

Figure 2.4: Energy spectra of the XXZ spin chain on 6 sites (reprinted from Ref. [13]).
The arrows indicate the points relating to the sl2 loop algebra symmetry. Color of thezzzz|: Magenta, |Stot| = 3; Green, |Stot| = 2; Red, |Stot| = 1;
lines indicates the value of |S tot
z
Blue, S tot = 0.

2.2

The spectral ﬂow of the XXZ spin chain on 6
sites.

Figure 2.4 is the energy spectral ﬂow for L = 6 , which is calculated numerically. Here,z| = 0 are, at
periodic boundary conditions are imposed: φ = 0. All the spectra for |S tot
least, double degenerated when L = 6.
The arrow in Fig. 2.4 indicates the points which is related to the sl 2 loop algebra
symmetry. As mentioned in Chapter 1, the sl 2 loop algebra symmetry leads to many
spectral degeneracies. The dimensions of degenerate eigenspaces are 12 at ∆ = ±0.5 and16 at ∆ = 0.

13

Chapter 3
Level Bifurcation
Branches such as bifurcations or multiple branches are observed in the ﬂow of eigenvalues
for the inhomogeneous transfer matrix of the XXZ spin chain [14]. The spectral ﬂow
for the XXZ anisotropic parameter ∆ is obtained by directly diagonalizing the transfer
matrix in the regime: −1 ≤ ∆ ≤ 1, where the inhomogeneous transfer matrix is real. The
branches are considered as novel spectral behaviors, which are not decomposed into simple
superpositions of level crossings. We may call them level bifurcations. We also observe thatno branches appear in the homogeneous case. The appearance of some of the branches are
conﬁrmed through the Bethe ansatz method. We note that an inhomogeneous transfer
matrix plays a central role in the exact calculation of thermodynamic quantities at ﬁnite
temperatures for the XXZ spin chain (for example, see Ref. [19]). We believe that thepresent study will shed some light on the spectral degeneracies in quantum systems.

3.1

Transfer matrix

Let us introduce the transfer matrix of the XXZ spin chain for both homogeneous andinhomogeneous cases [20, 21]. We deﬁne the L operator [21, 22] of the XXZ spin chain bysin(zIj + νσjz )sin(2ν)σj− 
Lj (z) = .
sin(2ν)σj+sin(zIj − νσjz )(3.1)

Here, I j is the identity matrix, σ j± is given by σj± = (σjx ± iσjy )/2, and the real parameterν is related to the parameter ∆ by

∆ = − cos(2ν)



π
0≤ν ≤
.
2

(3.2)

In this work, we consider only the massless regime of the XXZ spin chain: −1 ≤ ∆ ≤ 1.
All the nonzero matrix elements of L(z) are real and positive if z satisﬁes the condition:−ν < z < ν. We deﬁne the inhomogeneous transfer matrix τ (w; {ζ j }) acting on L sites
14

byτ (w; {ζj })= Tr(LL (w + ζL ) . . . L2 (w + ζ2 )L1 (w + ζ1 )).

(3.3)

Here, the parameters ζ1 , ζ2 , . . . , ζ L are called inhomogeneous parameters. The symbol {ζj }denotes the set of ζj ’s. When all ζj ’s take the same value, the transfer matrix is called
homogeneous. When they do not, it is called inhomogeneous. The Hamiltonian, Eq. (1.7),of the XXZ spin chain is given by the logarithmic derivative of the homogeneous transfer
matrix [20]:d∝ H + constant,(3.4)
ln τ 
dw
w=−ζ+ν
where the inhomogeneous parameters have the same value: ζj = ζ for j = 1, 2, . . . , L.

3.2

Eigenvalues of the transfer matrix obtained bynumerical diagonalization

Let us now discuss numerical results on the spectrum of the transfer matrix. Hereafter
we normalize the transfer matrix so that the eigenvalue is given by 1 for the “vacuum”
state in which all spins are up. We also set w = 0 for the transfer matrix. We consider
four cases for setting ζj ’s: three inhomogeneous cases and one homogeneous case. The
equations for the four cases, which correspond to Fig. 3.1, are the following:(a)ζj = kj ν(b)ζj = Lj · k0 ν

(0 < kj < 1; kj ’s are arbitrarily chosen.)(c) ζj = (−1)j · k0 ν(d)

ζj = k0 ν(homogeneous case)(3.5)

The real parts of eigenvalues of the transfer matrices τ (0; {ζj }) shown in Fig. 3.1
are obtained by numerically diagonalizing them for L = 6. In Fig. 3.1, we set k 0 = 0.8 in
Eq. (3.5). We observe branches in the inhomogeneous cases. On the other hand, there isno branch in the homogeneous case, as shown in Fig. 3.1(d).
Figure 3.2 shows that Case (b) in Eq. (3.5) is a special inhomogeneous case. In Case(b), two branching points lie on the same line, as shown in Fig. 3.2(b). On the contrary,
two branching points are not on the same line in Cases (a) or (c), as evident in Figs. 3.2(a)and (c).

15

2

(a)

Re Λ

1

0

-1

2/3

0.7

0.8

ν/(π/2)

0.9

1

2

(b)

Re Λ

1

0

-1

2/3

0.7

0.8

ν/(π/2)0.9

1Figure 3.1: The real part of the eigenvalues of the transfer matrices for L = 6 (partly
reprinted from Ref. [14]). Labels (a)-(d) correspond to the labels (a)-(d) in Eq. (3.5). The
spectral ﬂows in inhomogeneous cases have many branches, which do not appear in thehomogeneous case.

16

2

(c)

Re Λ

1

0

-1

0

0.2

0.4

0.6
ν/(π/2)

0.8

1

2

(d)

Re Λ

1

0

-1

2/3

0.7

0.8

ν/(π/2)

17

0.9

1

(a)

Re Λ

-0.4
-0.6
-0.8
-1
-1.2
0.7

0.8

ν/(π/2)

0.9

1

(b)

Re Λ

-0.4

-0.6
-0.8
-1
-1.2
0.7

0.8

ν/(π/2)0.9

1Figure 3.2: Graphs (a) and (b) show the magniﬁed spectral ﬂow in the sector of one
down spin, which has been given in Figs. 3.1(a) and 3.1(b), respectively (reprinted from
Ref. [14]). The solid lines denote single real eigenvalues, and the dashed lines denote realparts of complex eigenvalues.

18

Re Λ

-0.4
-0.6
-0.8
-1
-1.2
0.7

0.8

0.9

ν/(π/2)

1

Figure 3.3: The magniﬁed spectral ﬂow for L = 6 given by real solutions to the Bethe
ansatz equation, Eq. (3.6) (reprinted from Ref. [14]).

3.3

Eigenvalues of the transfer matrix obtained throughthe Bethe ansatz method

Let us conﬁrm the appearance of spectral branches by an independent method. First, weshow that the branches of the regular-interval case, i.e. Case (b) in Eq. (3.5), can also
be obtained by the Bethe ansatz method. For simplicity, we consider only the case of
one down spin. For the inhomogeneous case, according to Refs. [20, 21] the Bethe ansatz
equation with one down spin is given byLsin(t1 + ζj + ν)
j=1 sin(t1 + ζj − ν)

= 1.

(3.6)

The solution t1 leads to the eigenvalue of the transfer matrix:




L
L
sin(t1 − 2ν) sin(t1 + 2ν) sin(ζj + ν) +sin(ζj − ν) .
Λ(t1) =Anorm
sin(t1 ) j=1
sin(t1) j=1

1

(3.7)

Here, A norm denotes the normalization factor, which corresponds to the “vacuum” state
eigenvalue:
Anorm =L

j=1

sin(ζj + ν) +

L


j=1

sin(ζj − ν).

(3.8)

Now let us discuss numerical solutions to the Bethe ansatz equation, Eq. (3.6). We
introduce a function of rapidity t in the following:f (t) =L
j=1

sin(t + ζj + ν) −
19

L

j=1

sin(t + ζj − ν).

(3.9)

We can easily obtain real solutions to the equation: f(t) = 0, for instance, by the Newton
method. Substituting real solutions of Eq. (3.6) for Eq. (3.7), we obtain the eigenvaluesof the transfer matrix by the Bethe ansatz method. We note that Eq. (3.6) has L singular
points at the zeros of Lj=1 sin(t + ζj − ν) as a function of t.
We see in Fig. 3.3 that the same solid lines given in Fig. 3.2(b) are indeed obtainedwith the Bethe ansatz method. The numerical diﬀerences between the two methods areless than 10−7 for the eigenvalues of the transfer matrix. Thus, we may conclude thatthey are consistent. We have not reconstructed the dashed lines in Fig. 3.2(b) because
we can obtain only the real (not complex) solutions of Eq. (3.6) by the above method.

3.4

Discussion

Let us discuss the behavior of branching points for Case (b) in Eq. (3.5). Sending theinhomogeneous parameters continually to a ﬁxed value, we observe that the branching
points in the spectral ﬂow in Fig. 3.1(b) gradually move to the large-ν direction. Whenthe branching points approach the XXX point (ν = π/2), i.e. the homogeneous point, thebranches ﬁnally vanish as shown in Fig. 3.4. Figure 3.4 shows how the branching points
for Case (b) behave. Here, we rewrite Case (b) in Eq. (3.5) as




j−1
k1 k0 ν.
ζj = 1 −L(3.10)

When k1 = 1, Eq. (3.10) is eﬃciently equal to Case (b) in Eq. (3.5). When k 1 = 0,
Eq. (3.10) corresponds to the homogeneous case. In Fig. 3.4, the values of ν/(π/2) at thebranching points, which are calculated through the Bethe ansatz method, are shown as afunction of k1 . The solid lines and the broken lines in Fig. 3.4 correspond to the smallerν (ﬁrst) branching points and the larger-ν (second) branching points. Furthermore, inFig. 3.4(b) we show the ﬁtting lines to the values of ν/(π/2) at the branching points. The
ﬁtting function of the lines in Fig. 3.4(b) is given by
ν
1
,
=
π/2
1 + ak0 k1(3.11)

where a = 0.359 for the ﬁrst branching point, a = 0.178 for the second branching point,
and k0 = 0.8 for both points.
The spectral branches are important when we investigate the solutions of the Bethe
ansatz equations. In fact, it is nontrivial to solve the Bethe ansatz equations at a spectral
branch. Here, we recall Fig. 3.2(a). Two solid lines correspond to two real solutions of theBethe ansatz equation, Eq. (3.6), while the dashed line corresponds to the real part of twocomplex conjugate solutions of Eq. (3.6). Thus, two real solutions should be converted
into complex solutions at the branching point, where the Bethe ansatz equations may be
singular.
20

1

(a)

ν/(π/2)

0.9

0.8

0

k0=0.5, Firstk0=0.5, Second
k0=0.7, Firstk0=0.7, Second
k0=0.9, Firstk0=0.9, Second
0.2
0.4

0.6

0.8

1

k1
1

(b)

ν/(π/2)0.9

FirstSecond

0.8

0

0.2

0.4

0.6

0.8

1

k1

Figure 3.4: Branching points for Case (b) with k1 varying in Eq. (3.10). (a) Solid lines
are for the smaller-ν (ﬁrst) branching points; Broken lines are for the larger-ν (second)branching points. (b) Crosses are for the branching points when k0 = 0.8. The solid line
and the dashed line are the ﬁtting line to the crosses, where the ﬁtting function is given
by Eq. (3.11).

21

It is surprising that spectral branches appear in Case (c) in Eq. (3.5). The system has
the translational invariance in that case as well as the homogeneous case. We should saythat the expectation of Ref. [14] is incorrect: the appearance of spectral branches could be
related to the breakdown of the translational symmetry of the system. The expectation
has come from an ordinary interpretation. In other words, the level bifurcations are
nontrivial and extraordinary phenomena.

22

Chapter 4
Level Statistics
We have studied level statistics for XXZ spin chains with next-nearest-neighbor (NNN)couplings [15, 17] or with a random magnetic ﬁeld [16]. Evaluating the level-spacing
distribution, the spectral rigidity and the number variance, we conﬁrm the correspondence
between non-integrability and Wigner behavior in the spectrum. We show how the levelstatistics of the ﬁnite-size systems depends on the anisotropy parameter and the NNN
coupling parameter or the strength of the magnetic ﬁeld.

4.1

Preliminaries

Statistical properties of energy levels have been studied for various physical systems interms of the random matrix theory (RMT). In the RMT, we consider a large Hamiltonian whose elements are random variables with given probability laws. The statisticalbehavior of the energy levels depends only on the symmetry of the Hamiltonian [23].
For systems with time-reversal invariance, the Gaussian orthogonal ensemble (GOE) isappropriate for the description of the Hamiltonian matrices. For systems having invariance under time-reversal but no rotational symmetry, the Gaussian symplectic ensemble(GSE) is appropriate. For systems without time-reversal invariance, the Gaussian unitaryensemble (GUE) is appropriate. The RMT analysis was performed also in the view ofsemi-classical pictures of quantum many-body systems [24, 25]. When the classical motion is integrable, level-clustering can occur. Then the level-spacing distribution is thenegative exponential characteristic of Poisson processes [24]. Chaotic motion in the classical system is associated with level repulsion and with the level statistics characteristicof GOE [25].
For quantum many-body systems, the RMT analysis has been applied to characterizing quantum chaos and to investigating the integrability of a system. For quantum
spin systems, we adopt a deﬁnition of integrability by the Bethe ansatz: an integrable
model is exactly solvable by the Bethe ansatz. Since a pioneering work [26], the following
23

conjecture has been widely accepted: If a given Hamiltonian is integrable by the Bethe
ansatz, the level-spacing distribution should be described by the Poisson distribution:
PPoi (s) = exp(−s).

(4.1)

If it is nonintegrable, the level-spacing distribution should be given by the Wigner distribution:


πs
πs2
PWig (s) =
.
(4.2)
exp −
2
4
In fact, the above conjecture has been numerically conﬁrmed for many quantum spin
systems such as correlated spin systems [15, 17, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35] anddisordered spin systems [16, 36, 37, 38, 39]. In the Anderson model of disordered systems,
PPoi (s) and PWig (s) characterize the localized and the metallic phases, respectively [40].
Here we note that the Wigner distribution, Eq. (4.2), is commonly known as the
Wigner surmise for GOE. In this chapter, we use the word “Wigner behavior” instead of“GOE behavior”.
It should be important to study statistical properties of energy levels for XXZ spin
chains. The XXZ spin chains are closely related to the most important families of integrable spin systems through the integrable point as mentioned in Chapter 1. On theother hand, XXZ spin chains are also extended into nonintegrable systems quite naturallywith NNN couplings or a random magnetic ﬁeld. We can therefore expect that the RMTanalysis of the XXZ spin chains could be important in discussing level statistics for otherquantum systems that would have some connection to an integrable system.
The organization of this chapter is as follows: In Sec. 4.2, we describe our models andnumerical methods. We explain the details especially about desymmetrization and unfolding. Deﬁnitions of level-spacing distribution, spectral rigidity and numerical varianceare also described in the section. In Sec. 4.3, we discuss level statistics of integrable XXZ
spin chains in both a standard case and a special symmetrical case. In Sec. 4.4, the characteristic behavior of level statistics is discussed for NNN coupled XXZ spin chains. First
we conﬁrm the correspondence between non-integrability and Wigner behavior. Then weshow how the behavior of level statistics changes due to mixed symmetry and ﬁnite-size
eﬀects. We also show that the characteristic behavior of level statistics does not dependin the energy range. In Sec. 4.5, the behavior of level statistics is discussed for XXZ
spin chains under random magnetic ﬁeld. We show how the behavior of level statistics
depends on the anisotropy parameter, the magnetic ﬁled and the system size. We also
discuss the correspondence between the Anderson model of disordered systems and ourmodel. Finally, we summarize this chapter in Sec. 4.6.

24

(a)

=
J1

J2

J1

J2

h5

(b)

h4

5

4
3 h3

N
hN

2

1

h2

h1

Figure 4.1: Schematic pictures of the models. (a) XXZ spin chain with next-nearestneighbor couplings. It is also called “zig-zag chain”. (b) XXZ spin chain with a randommagnetic ﬁeld.

4.2

Models and numerical procedure

4.2.1

ModelsLet us give the Hamiltonian of a XXZ spin chain with NNN couplings on L sites by
Hnnn = J1

L

j=1

y
x
z
(SjxSj+1+Sjy Sj+1+∆1Sjz Sj+1
)+J2

L


y
x
z
(Sjx Sj+2+Sjy Sj+2+∆2 Sjz Sj+2). (4.3)

j=1

We also give that of a XXZ spin chain under random magnetic ﬁeld byHrandom = J

L 

j=1



y
x
z
Sjx Sj+1+ Sjy Sj+1
+ ∆Sjz Sj+1
+

L

j=1

hj Sjz .

(4.4)

Here, S α = 12 σ α and (σ x, σ y , σ z ) are the Pauli matrices; h j is random magnetic ﬁeld along
the z axis at site j; the periodic boundary conditions are imposed. The random magnetic
ﬁeld hj ’s are uncorrelated random numbers with a Gaussian distribution: hj  = 0 andhn hm  = h2 δnm. The schematic pictures of the above two models are shown in Fig. 4.1.
The Hamiltonian (4.3) is nonintegrable when the NNN coupling parameter J 2 isnonzero, while it is integrable when J 2 vanishes. Here we note that Eq. (4.3) is coincides
with the NNN coupling Heisenberg chain when ∆ 1 = ∆2 = 1. We should also note thatthe Hamiltonian (4.4) is integrable when h = 0, while it is not when h = 0.

25

4.2.2Desymmetrization

Let us discuss desymmetrization of the Hamiltonians of the XXZ spin chains. Whenperforming calculation on level statistics, one has to separate the Hamiltonian matrices
into some sectors; in each sector, the eigenstates have the same quantum numbers. The
NNN coupled XXZ chains are invariant under spin rotation around the z axis, translation,
reﬂection, and spin reversal. Therefore we consider quantum numbers for the total S zz(Stot), the total momentum Ktot , the parity, and the spin reversal. However the total
momentum Ktot is invariant under reﬂection only when Ktot = 0 or π. Thus the desymz
ismetrization according to parity is performed only when Ktot = 0 or π. Similarly, S totzinvariant under spin reversal only when Stot = 0. Thus the desymmetrization accordingz= 0. In the case of a XXZ spin chain underto spin reversal is performed only when Stotz
as a quantum number.
random magnetic ﬁeld, we have to consider only Stot
It is convenient to use a momentum-based form for the Hamiltonian when we discuss
the symmetries of the model as well as when we calculate the eigenvalues. To obtain theform, we perform the Jordan-Wigner and the Fourier transformations on the original
Hamiltonian. Some details are explained in Appendix B.
We now discuss some details of the spin reversal operation. The spin reversal operation on the spin variable of the jth site is deﬁned by
Sj± → Sj∓ ,Sjz → −Sjz .

(4.5)

Here, S j± = (Sjx ± iSjy )/2. Let M denote the number of down-spins in a given sector. Thezz
is given by Stot
= L/2 − M . We can show that invalue of the total spin operator Stotz
= L/2 − M , the transformation of Eq. (4.5) corresponds to the following
a sector of Stot
operation on momentum-based fermion operators:ĉ†k → −ĉπ−k ,ĉk → −ĉ†π−k .

(4.6)

Here, ĉ†k and ĉk are the creation and annihilation operators of free fermions with momentum k. We note that the value of momentum k is given by (2π/L)× (an integer) for odd
M and (2π/L)× (a half-integer) for even M . An explicit derivation of the transformation
Eq. (4.6) from the deﬁnition Eq. (4.5) is given in Appendix C.
For desymmetrizing the Hamiltonian with respect to spin reversal symmetry, it
is useful to know how the vacuum state transforms under the spin reversal operation
expressed in terms of the momentum-based fermion operators. Let us denote by |0
the vacuum state where there is no down-spin. Under the spin reversal operation, it
transforms up to a phase factor AL as follows|0 → AL ĉ†q1 ĉ†q2 · · · ĉ†qL |0(4.7)

Here q j denotes momentum (2π/L)j for j = 1, 2, . . . , L, when M is odd, and (2π/L)(j −1/2) for j = 1, 2, . . . , L, when M is even. We can show that the phase factor AL is given
26

by
AL =

1



LL/2 P ∈SLP exp(−iL
j=1

jkP j ) .

(4.8)

Here SL denotes the set of permutations on L elements, P the sign of permutation P .
Furthermore, we can calculate the phase factor A L as follows



AL =

(−1) forM = 2
+1
 (−1)
for M = 2 + 1(4.9)

The derivation is given in Appendix D.
It is sometimes convenient to use the following in stead of Eq. (4.6):ĉ†k → ĉπ−k ,ĉk → c†π−k .

(4.10)The Hamiltonians are invariant not only for Eq. (4.6) but also for Eq. (4.10). The form
Eq. (4.10) has an advantage that we do not need to consider the phase factor −1 thatappears in Eq. (4.6).
Making use of the spin reversal operation expressed in terms of the fermion basiszEq. (4.10), we have desymmetrized the Hamiltonian matrix in the sector S tot= 0 withzrespect to spin reversal symmetry. For a given vector with S tot = 0, we calculate how
it transforms under the operation Eq. (4.10). If it is not a singlet and transforms into a
diﬀerent vector, then we combine the pair into an eigenvector of the operation Eq. (4.10).

4.2.3

Unfolding

To ﬁnd universal statistical properties of the Hamiltonians, one has to deal with unfolded
eigenvalues instead of raw eigenvalues. The unfolded eigenvalues are renormalized values,
whose local density of states is equal to unity everywhere in the spectrum. In this work,
the unfolded eigenvalues x i’s are obtained from the raw eigenvalues Ei in the following
method. Let us deﬁne the integrated density of states asn(E) =N
i=1

θ(E − Ei ).

(4.11)

Here, θ(E) is the step function and N is the number of the eigenvalues. We choose somepoints of coordinates: (Ei , n(Ei )) for i = 1, 21, 41, · · · , N . The average of the integrated
density of states n(E) is approximated by the spline interpolation through the chosenpoints. The unfolded eigenvalues are deﬁned as
xi = n(Ei ).

(4.12)In Fig. 4.2, the schematic picture of the integrated density of states shows how to obtain
unfolded eigenvalues from raw eigenvalues.
27Integrated density of states

9 x9
8
x8
7 x7
6 x6
5 x5
4 x4
3 x3
2 x2
1 x1
0
E1

n(E)
n(E)

E2 E3 E4 E5
E6 E7 E8 E9EnergyFigure 4.2: How to obtain unfolded eigenvalues xj ’s from raw eigenvalues Ej ’s. The
integrated density of states n(E) is deﬁned by Eq. (4.11). The average of the integrated
density of states n(E) here is just an example.

4.2.4

Level-spacing distribution, spectral rigidity and number
variance

To analyze spectral properties, in this chapter, we calculate three quantities: level-spacing
distribution P (s), spectral rigidity ∆ 3(l), and number variance Σ2 (l). The level-spacing
distribution is the probability function P (s) of nearest-neighbor level-spacing s = x i+1 −xi,
where xi ’s are unfolded eigenvalues. The level-spacing distribution is calculated over thewhole spectrum of unfolded eigenvalues unless we specify the range. The spectral rigidity
is given by ε0 +l/2
12[Nu (ε) − aε − b] dε,(4.13)min∆3 (l) =l a,b ε0 −l/2ε0where Nu (ε) = Σi θ(ε − εi ) is the integrated density of unfolded eigenvalues and ε0denotes an average over ε0 . The average is done on the whole spectrum except about 15
levels on each side. The expression of ∆ 3(l) gives the least square deviation of Nu (ε) fromthe best ﬁt straight line in an interval l. The number variance is given by
2Σ (l) =



Nu







l
l
ε0 +
− Nu ε0 −
−l
2
2

2 

,

(4.14)

ε0

where ε0 denotes an average over ε0 [34]. The average is done on the whole spectrum
except about 10 levels on each side.
For NNN coupled chains, we calculate the spectra for L = 18. The matrix size is28

Table 4.1: Matrix size of the subspaces where calculation is performed for NNN coupled
XXZ spin chains (L = 18) in this work.
z
Stot
and Ktot
Matrix sizez
Stot = 0 and Ktot = 2π/181364 × 1364
z
Stot = 1 and Ktot = 01282 × 1282
z
Stot = 0 and Ktot = 0 (not considering the spin reversal) 1387 × 1387

Table 4.2: Matrix size of the largest subspace and the number of samples calculated inthis work for XXZ spin chains under random magnetic ﬁeld.
Lattice size
8
10
12
14

Matrix size Number of samples
70 × 70
10000
252 × 2523000
924 × 924
1000
3432 × 3432
300

zgiven by the following (See also Table 4.1): 1364 × 1364 for Stot
= 0 and Ktot = 2π/L,zwhere L is the number of sites; 1282 × 1282 for Stot = 1 and Ktot = 0; 1387 × 1387z
= 0 and Ktot = 0 when desymmetrization is performed except for spin reversal.
for Stot
For XXZ spin chains under a random magnetic ﬁeld, we consider the largest subspacez= 0. The largest sectors for the lattice size L =8, 10, 12, 14 have 70, 252, 924,
Stot
3432 eigenvalues, respectively (see Table 4.2). We have calculated 10000, 3000, 1000, 300
samples of P (s) and ∆3 (l) for L =8, 10, 12, 14, respectively (see Table 4.2), and averaged
the samples for each L. To calculate the eigenvalues, we use standard numerical methods,
which are contained in the LAPACK library.
Numerical calculations are performed for the XY -like region, |∆| < 1, where ∆is the anisotropic parameter. It may be interesting to study for the Ising-like region,|∆| > 1, because there exist Ising-like magnets. For example, CsCoBr 3 and CsCoCl3
are quasi-1D Ising-like antiferromagnets with ∆ ∼ 10 [10]. For ∆1, however, levelstatistics is not reliable because energy spectra have some large gaps relative to ∆ andthe above unfolding method is invalid.

29

1
0.8

P(s)

0.6
0.4
0.2
0
0

0.5

1

1.5

s

2

2.5

3

0.5
0.4

∆3(l)

0.3
0.2
0.1
0
0

5

l

10

15

20

1.4
1.2

Σ2(l)

1
0.8
0.6
0.4
0.2
0
0

2

4

l

6

8

10

Figure 4.3: Level-spacing distribution, spectral rigidity and numerical variance of thez= 1, Ktot = 0
integrable XXZ spin chain (J 2 = 0) for L = 18, J1 = 1, ∆1 = 0.8, Stot
(reprinted from Ref. [17]). Broken line shows the Poisson distribution.

30

4.3

Integrable XXZ spin chain

Let us conﬁrm the Poissonian behavior for the generic case of the integrable XXZ spin
chain. When J2 = 0, the level-spacing distribution P (s) mostly show the Poisson distribution as shown in Fig. 4.3. The spectral rigidity ∆3 (l) and the number variance Σ2 (l)also shows Poissonian behavior. We conﬁrmed numerically the standard result for some
generic values of the XXZ coupling ∆ 1 .
On the other hand, extraordinary behavior can be seen for the special case of theintegrable XXZ spin chain that has the sl 2 loop algebra symmetry. Here we recall that
the XXZ spin chain is integrable when the NNN coupling J 2 vanishes, and also that thesl2 loop algebra symmetry exists when q is a root of unity. Here the anisotropy ∆ 1 isrelated to q through ∆1 = (q + q −1)/2. For instance, when the parameter q is given byexp(iπ/3), we have ∆1 = 0.5.
The level spacing distribution P (s) and the spectral rigidity ∆ 3(l) of the integrable
XXZ spin chain are shown for ∆1 = 0.5 in Fig. 4.4. In Fig. 4.4 (a) we do not perform desymmetrization with respect to spin reversal, while in Fig. 4.4 (b) we plot thelevel-spacing distribution P (s) and the spectral rigidity ∆ 3(l) after we desymmetrize the
Hamiltonian with respect to the spin reversal operation.
We observe that there still remain many degeneracies associated with the sl 2 loop algebra symmetry, even after desymmetrizing the Hamiltonian with respect to spin reversal
symmetry. The level-spacing distribution P (s) has a small peak at s = 0 in Fig. 4.4 (b).
Furthermore, the slopes of ∆3(l) shown in the insets of Figs. 4.4 (a) and 4.4 (b) are larger
than that of Poissonian behavior. However, the numerical result does not necessarily givea counterexample to the conjecture of RMT. The level statistics might show Poissonian
behavior, if we completely desymmetrize the Hamiltonian matrix in terms of the sl 2 loop
algebra symmetry.

4.4

Next-nearest-neighbor coupled XXZ spin chains

We now numerically discuss the characteristic behavior of level statistics for next-nearestneighbor (NNN) coupled XXZ spin chains in this section.

4.4.1

z
Wigner behavior for Stot=0z
For a sector of Stot
= 0, we numerically discuss the characteristic behavior of level statistics
on the XXZ spin chains. We should remember that the system is invariant under spinz
= 0. Here we note that spin reversal symmetry has not been
reversal only when Stotconsidered explicitly in previous studies of level statistics for various quantum spin chains.
In some sense, desymmetrizing the Hamiltonian with respect to spin reversal symmetry

31

(a)

∆3(l)

4

P(s)

3

0.4
0.3
0.2
0.1
0

5

10

15

2

2.5

10

15

2

2.5

l

2

1

4

0.5

(b)

P(s)

3

1

∆3(l)

0
0

s

1.5

0.4
0.3
0.2
0.1
0

5

3

l

2

1

0
0

0.5

1

s

1.5

3

Figure 4.4: Level-spacing distribution of the integrable XXZ spin chain (J 2 = 0) forz= 0, Ktot = 2π/L. (a) Desymmetrization with respect
L = 18, J1 = 1, ∆1 = 0.5, Stotto spin reversal is not performed. (b) Desymmetrization with respect to spin reversal is
performed. The inset is spectral rigidity for each case (reprinted from Ref. [17]). Brokenlines corresponds to Poissonian behavior, and solid lines Wigner behavior.

32

has been avoided due to some technical diﬃculty. For XXZ spin chains, Level statisticsz= 0, where it is not necessary to consider spinhas been discussed only for sectors of Stotreversal symmetry.
z= 0 if we
Let us show explicitly such a case that Wigner behavior appears for Stotconsider spin reversal symmetry. In Fig. 4.5 we have obtained the numerical results forlevel statistics such as the level-spacing distribution P (s), the spectral rigidity ∆ 3 (l) andz= 0. Here we note thatthe number variance Σ2 (l), for the sector of Ktot = 2π/L and Stotin the sector the parity invariance does not exist and we focus on spin reversal symmetry.
The numerical results of level statistics shown in Fig. 4.5 clearly suggest Wigner behavior.
The curve of the Wigner distribution ﬁts well to the data of the level-spacing distribution
P (s). The plots of the spectral rigidity ∆ 3 (l) are consistent with the curve of Wigner
behavior as shown in Fig. 4.5. It is also the case with the number variance Σ2 (l).

4.4.2

How does non-Wigner behavior appear?

Unexpected non-Wigner behavior has been reported in Ref. [15] for level-spacing distributions of the NNN coupled XXZ chains. Let us discuss the reason why it was observed,
considering both mixed symmetry and ﬁnite-size eﬀects. There are two types of nonWigner proﬁles for the nonintegrable systems: one is given by almost the numerical
average of the Poisson and the Wigner distributions, and another one is rather close to
the Poisson distribution. The proﬁles of the ﬁrst type appear in various cases [15], suchas the case of ∆2 = 0.5. Similar non-Wigner behavior has been observed for a circularbilliard when the angular momentum Lz = 0, and for an interacting two-electron systemwith the Coulomb interaction in a quantum billiard when L z = 0 [41]. The proﬁles of thesecond type appear in particular for the case of ∆1  1 or ∆2  1. We may call the latter
Poisson-like behavior rather than simple non-Wigner behavior. Both types of non-Wignerz= 0, which is the largest sector
distributions have been observed in the subspace of Stotof the Hamiltonian matrix of Eq. (4.3). Here we note that the observations in Ref. [15]zz
= 0 do not contradict with the Wigner behavior of Ref. [31] for Stot= 0, where
for Stot
the level-spacing distributions of similar XXZ chains have been discussed.
z= 0 or 1 andWe show level-spacing distributions in Fig. 4.6 for the four cases: Stotzshould be∆1 = ∆2 = 0.5 or 0.98. The numerical results suggest that the value of Stot
important as well as the anisotropy parameters, ∆1 and ∆2, in the observed non-Wigner
behavior of the level-spacing distributions. When ∆ 1 = ∆2 = 0.5, Wigner behaviorzz= 1, while the non-Wigner behavior was observed for Stot
= 0. We haveappears for Stotzalso checked that Wigner behavior appears for Stot = 2. Furthermore, we have conﬁrmedz-dependence of the level-spacing distribution is valid for some values of K tot .
that such Stotz
, Ktot and the parity when
Here we have desymmetrized the Hamiltonian according to S tot
it exists, but not to the spin reversal. Here we note that the parity invariance exists only
33

1
0.8

P(s)

0.6
0.4
0.2
0
0

0.5

1

s

1.5

2

2.5

3

0.5

∆3(l)

0.4
0.3
0.2
0.1
0
0

5

l

10

15

20

1.4
1.2

Σ2(l)

1
0.8
0.6
0.4
0.2
0
0

2

4

l

6

8

10

Figure 4.5: Level-spacing distribution P (s), spectral rigidity ∆ 3(l), and number variancez= 0,
Σ2 (l) of the NNN coupled chain for L = 18, J 2 /J1 = 0.5, ∆1 = ∆2 = 0.5, Stot
Ktot = 2π/L under complete desymmetrization (reprinted from Ref. [17]). Broken lines
correspond to Poissonian behavior, and solid curves Wigner behavior.

34

1

1
z
tot

z

S =0
∆1=∆2=0.5

0.8

Stot =1∆1=∆2=0.5

0.8

P(s)

0.6

P(s)

0.6
0.4

0.4

0.2

0.2

0
0

0.5

1

s

1.5

2

2.5

0
0

3

1

0.5

1

s

1.5

2

2.5

3

1
z
tot

S =0
∆1=∆2=0.98

0.8

Sztot =1∆1=∆2=0.98

0.8

P(s)

0.6

P(s)

0.6
0.4

0.4

0.2

0.2

0
0

0.5

1

s

1.5

2

2.5

3

0
0

0.5

1

s1.5

2

2.5

3

Figure 4.6: Level-Spacing distribution P (s) of the NNN coupled chain for L = 18, J 2 /J1 =1, Ktot = 0 (reprinted from Ref. [17]). Here, desymmetrization is performed except for
spin reversal. Broken lines, the Poisson distribution; solid curves, the Wigner distribution.
In each of the four distributions P (s) vanishes at s = 0: there is no degeneracy among
energy levels.

35

0.5

∆3(l)

0.4
0.3
0.2

Stotz =0, ∆1=∆2=0.98
Stotz =1, ∆1=∆2=0.98
Stotz =0, ∆1=∆2=0.5
Stot=1, ∆1=∆2=0.5
z

0.1
0
0

5

l

10

15

20

1.4
1.2

Σ2(l)

1
0.8
0.6
0.4
0.2
0
0

2

4

l

6

8

10

Figure 4.7: Spectral rigidity ∆ 3(l) and number variance Σ2 (l) of the NNN coupled chain
for L = 18, J2 /J1 = 1, Ktot = 0 (reprinted from Ref. [17]). Here, desymmetrization is
performed except for spin reversal. Broken lines, Poissonian behavior; Solid lines, Wigner
behavior.

for sectors with Ktot = 0 or π when L is even.
z= 0 and ∆1 = ∆2 = 0.5 shown in
The non-Wigner behavior observed for the case Stot
Fig. 4.6 should be due to mixed symmetry. We have not performed the desymmetrizationz= 0. Let us recall the Wigner
according to spin reversal which is important only for Stotzbehavior for Stot = 0 shown in Fig. 4.5. It is for Ktot = 2π/L. However, Wigner behaviorz= 0 where we have to perform theshould appear also for the sector of Ktot = 0 and Stotdesymmetrization according to both parity and spin reversal simultaneously. Because thebehavior of level statistics should be independent of K tot .
The Poisson-like behavior for the case ∆1 = ∆2 = 0.98 should be dominated byz= 1 of Fig. 4.6, the Poisson-like behavior appears
the ﬁnite-size eﬀects. In fact, for S tot
when ∆1 = ∆2 = 0.98, while Wigner behavior appears when ∆1 = ∆2 = 0.5. We have
conﬁrmed that such tendency does not depend on the value of Ktot : we see it not only
for Ktot = 0 but also for Ktot = 0. We discuss the ﬁnite-size eﬀects on this model later.
The observations of the level-spacing distributions can also be conﬁrmed by investigating spectral rigidity ∆ 3(l) and number variance Σ2 (l). In Fig. 4.7, ∆ 3(l) and36

J2/J1
η= 0
η=0.2
η=0.4
η=0.6
η=0.8
η= 1

0.9

Wigner

0.7
(a)

0.5

Poisson

0.3
(b)
0.1
-1

-0.8

-0.6
-0.4
∆ (=∆1=∆2)

-0.2

0

0

0.2

0.4
0.6
∆ (=∆1=∆2)

0.8

1Figure 4.8: The diagram of contour lines of η for the NNN coupled chain with L = 18,z= 1 and Ktot = 0 (reprinted from Ref. [17]). Roughly speaking,in the sector of Stot
the area above the red long dashed line, Wigner behavior; the area below the blue solidline, Poissonian behavior. The points (a) and (b) correspond to Figs 4.9(a) and 4.9(b),respectively.

z= 1 andΣ2 (l) are shown for the four cases corresponding to those of Fig. 4.6. For Stotz
∆1 = ∆2 = 0.5, Wigner behavior appears. For Stot = 0 and ∆1 = ∆2 = 0.5, an intermediate behavior appears, which is close to the average between Wigner and Poissonianzz= 0 and Stot
= 1, Poisson-like behavior
behaviors. For ∆1 = ∆2 = 0.98, both for Stot2
appears for ∆3(l) and Σ (l).

4.4.3

Finite-size eﬀects on the level spacing distribution

Let us explicitly discuss ﬁnite-size eﬀects appearing in level statistics. Hereafter we discussz= 1 so that we do not have to perform the desymmetrizationlevel statistics for the case S tot
according to spin reversal. They should be important in the Poisson-like or non-Wigner
behavior observed in level statistics for the completely desymmetrized XXZ Hamiltonians.
There are two regions in which ﬁnite-size eﬀects are prominent: A region where J 2 isclose to zero and another region where ∆1 and ∆2 are close to 1. In the former region,
quantum integrability appears through ﬁnite-size eﬀects, and the characteristic behaviorof level statistics becomes close to Poisson-like behavior. In the latter region, Poisson-likebehavior appears due to the symmetry enhancement at the point of ∆1 = ∆2 = 1, where
the U (1) symmetry of the XXZ spin chain expands into the spin SU (2) symmetry.
Let us now discuss how the degree of non-Wigner behavior depends on the anisotropy
parameters, ∆1 and ∆2, and the NNN coupling, J 2 . For simplicity we set ∆ 1 = ∆2 anddenote it by ∆, and we also consider the ratio of J2 /J1 . We express the degree of non-

37

Wigner behavior by the following parameter: s0[P (s) − PWig (s)]ds,
0 [PPoi (s) − PWig (s)]dsη =  s00

(4.15)

where s0 = 0.4729 · · · is the intersection point of PPoi (s) and PWig (s) [16, 37]. We have
η = 0 when P (s) coincides with PWig (s), and η = 1 when P (s) coincides with PPoi (s).
The diagram of contour lines of η is shown in Fig. 4.8. We have calculated them for thearea −0.98 ≤ ∆ ≤ −0.02, 0.02 ≤ ∆ ≤ 0.98, and 0.02 ≤ J2/J1 ≤ 1, where ∆ = ∆1 = ∆2.
The contour lines of η show that behavior close to Wigner one appears in a large
region, while the Poisson-like behavior appears in a narrow region along the line of J2 /J1 =0 and that of ∆1 = ∆2 = 1. The Poisson-like behavior is dominated by ﬁnite-size eﬀects
and hence should vanish when L → ∞. This expectation should be consistent with thesuggestion in Ref. [31] that an inﬁnitesimal integrability-breaking term (the NNN term ofEq. (4.3) in this paper) would lead to Wigner behavior. Here we remark that the phase
diagrams of the ground state [42, 43] are totally diﬀerent from the diagram of contour
lines of η. It is due to the fact that level statistics reﬂects highly excited states rather
than the ground state.
When ∆1 = ∆2 = 0.98, Poisson-like behavior appears in level statistics due to someﬁnite-size eﬀects. It will be explicitly shown in Figs. 4.6 and 4.7. When ∆ 1 = ∆2 = 1,
Eq. (4.3) coincides with the Heisenberg chain, which has the spin SU (2) symmetry. Some
degenerate energy levels at ∆ 1 = ∆2 = 1 can become nondegenerate when ∆1 and ∆2 arenot equal to 1. The diﬀerence among the nondegenerate energy levels should be smaller
than the typical level spacing when ∆ 1 and ∆2 are close to 1. The typical level spacing,which is of the order of 1/L, should become large when the system size L is small. Thus,
the Poisson-like behavior should practically appear in level statistics. We note that for
the Heisenberg chain Wigner behavior appears in the level-spacing distribution when wedesymmetrize the Hamiltonian with respect to the spin SU (2) symmetry [27, 28].

4.4.4

Homogeneity of the characteristic behavior of level statistics throughout the spectrum

Let us discuss that for the XXZ spin chains the characteristic behavior of level statistics
does not depend on the energy range of the spectrum. In Fig. 4.9, we show level-spacing
distributions evaluated at points (a) and (b) shown in the diagram of Fig. 4.8. They areevaluated for three diﬀerent energy ranges. The distributions shown in Fig. 4.9(a) give
Wigner behavior, while the distributions of Fig. 4.9(b) are close to Poisson-like behavior.
Let us explain the three diﬀerent energy ranges shown in Fig. 4.9. Red histograms
show the level-spacing distributions evaluated for all levels, while green bars show thoseevaluated only for the 1/3 of all levels around the center, and blue crosses for the 10% of
all levels located from each of the two spectral edges.
38

1

(a)
0.8

P(s)

0.6
0.4
0.2
0

0

0.5

1

s

1.5

2

2.5

3

1

(b)
0.8

P(s)

0.6
0.4
0.2
0

0

0.5

1

s

1.5

2

2.5

3

Figure 4.9: Level-spacing distribution of the NNN coupled chain for L = 18, ∆ 1 = ∆2 =z
= 1, and (a) J2/J1 = 0.5; (b) J2 /J1 = 0.2 (reprinted from Ref. [17]). Red
0.5, Stothistograms are for all levels; green bars, 1/3 of all levels around the center; blue crosses,
10% of all levels from each of the two edges. Solid and broken lines show the Wigner andPoisson distributions, respectively.

39

Quite interestingly, the distributions evaluated for the diﬀerent energy ranges arequite similar to each other. Here we should note that some statistical ’noises’ exist but they
are rather small. This makes a remarkable contrast between the level-spacing distributionsof the XXZ spin chains and that of the Anderson model of disordered systems. Here werecall that even in the metallic phase the level-spacing distribution of the Anderson modelshows Poissonian behavior if we evaluate it around the two edge regions of the energy
spectrum.

4.5

XXZ spin chains under random magnetic ﬁeld

Now we discuss the characteristic behavior of the XXZ spin chains under random magnetic
ﬁeld. Under random magnetic ﬁeld, the system is nonintegrable.

4.5.1

Dependence on the anisotropy parameter ∆Depending on the anisotropic parameter ∆, the level-spacing distribution P (s) and thespectral rigidity ∆ 3(l) change between the Wigner behavior and the Poissonian behavior
as shown in Fig. 4.10, where L = 14 and h/J = 0.5. When ∆ = 0, P (s) almost coincides
with PPoi (s) although the system is nonintegrable due to the random magnetic ﬁeld. As∆ increases, P (s) rapidly changes to PWig (s). The behavior of ∆3(l) also changes fromPoissonian to Wigner behavior as ∆ increases.
Let us explain the Poissonian behavior of ∆ = 0 in terms of the Anderson localization. The Hamiltonian (4.4) can be mapped into a model of interacting one-dimensional(1D) free fermions under random potential (see Appendix B):




 †
J L−1
2(cj cj+1 + c†j+1 cj ) − (−1)M (c†L c1 + c†1 cL )
H =
4 j=1

+ ∆
+

L
J
(4c†j cj c†j+1 cj+1 − 4c†j cj + 1)
4 j=1

L

j=1



hj



1
− c†j cj .
2

(4.16)

Here, L is the number of sites; M is the number of fermions; c †j and cj are the creation
and annihilation operators of fermions on the jth site, respectively. And the Andersonmodel of noninteracting disordered fermions is given byH=jεj c†j cj +


i,j

V (c†i cj + c†j ci ),(4.17)

where εj is the random potential at the jth site; V is a constant hopping integral; i, j
denotes summation over nearest-neighbor sites. One can ﬁnd that Eq. (4.16) for ∆ = 0
40

1

∆=1
∆=0.5
∆=0.2
∆=0.1
∆=0

0.8

P(s)

0.6
0.4
0.2
0
0

0.5

1

s

1.5

2

2.5

3

0.5

∆3(l)

0.4
0.3
0.2
0.1
0
0

5

l

10

15

20

Figure 4.10: Level spacing distribution P (s) and spectral rigidity ∆ 3 (l) of the XXZ chains
under random magnetic ﬁeld for L = 14, h/J = 0.5 (partly reprinted from Ref. [16]).
Broken lines, Poissonian behavior; solid lines, Wigner behavior.

41

1
1.2

h/J=0.1
h/J=0.5
h/J=1
h/J=2

1

h/J=0.1
h/J=0.5
h/J=1
h/J=2

0.8

0.8

P(s)

0.6

L=8

P(s)

0.6

L=10

0.4

0.4
0.2

0.2
0
0

0.5

1

s

1.5

2

2.5

0
0

3

1

1

s

1.5

2

2.5

3

1

h/J=0.1
h/J=0.5
h/J=1
h/J=2

0.8

h/J=0.1
h/J=0.5
h/J=1
h/J=2

0.8

P(s)

0.6

P(s)

0.6

L=12

L=14

0.4

0.4

0.2

0.2

0
0

0.5

0.5

1

s

1.5

2

2.5

3

0
0

0.5

1

s

1.5

2

2.5

3

Figure 4.11: Level-spacing distributions of the XXZ chains under random magnetic ﬁeld
for L = 8, 10, 12, 14; h/J =0.1, 0.5, 1, 2; ∆ = 0.5 (reprinted from Ref. [16]).

corresponds to Eq.(4.17) of the 1D case. It is known that localization always occurs in
the 1D case, while the 3D Anderson model has the metallic phase and the localized phase.
Here we recall that the metallic phase corresponds to P Wig (s) and the localized phase to
PPoi (s). Thus, the observed Poissonian behavior for ∆ = 0 is consistent with the Anderson
localization.
The Hamiltonian (4.4), namely Eq. (4.16), for ∆ = 0 corresponds to the interacting
1D fermions under random potential. Then, does the Wigner-like behavior of P (s) for∆ = 0 suggest that the interaction among fermions should break the Anderson localizationin 1D chains? The answer should be “No”. According to Ref. [44], the ground state isin the localized phase at least when 0 < ∆ < 1. The characteristic behavior of levelstatistics is inﬂuenced mainly by highly excited states, not by the ground states.

4.5.2

Dependence on the magnetic ﬁeld

We now discuss how the level-spacing distribution P (s) depends on the random magnetic
ﬁeld h. We consider only the case of ∆ = 0. In Fig. 4.11, the graphs of P (s) are shown
42

0.5

L=8

0.4

0.3

∆3(l)

∆3(l)

0.4

0.5

0.2

h/J=2
h/J=1
h/J=0.5h/J=0.1

0.1
0
0

5

l

10

15

h/J=2
h/J=1
h/J=0.5h/J=0.1
5

l

10

15

20

0.5

L=12

0.4

∆3(l)

∆3(l)

0.2

0
0

20

0.3
0.2

h/J=2
h/J=1
h/J=0.5h/J=0.1

0.1
0
0

0.3

0.1

0.5
0.4

L=10

5

l

10

15

L=14

0.3
0.2

h/J=2
h/J=1
h/J=0.5h/J=0.1

0.1
20

0
0

5

l

10

15

20

Figure 4.12: Spectral rigidity ∆ 3(l) of the XXZ chains under random magnetic ﬁeld for
L = 8, 10, 12, 14; h/J =0.1, 0.5, 1, 2; ∆ = 0.5.

for some values of h/J and L, where ∆ = 0.5. We ﬁrst consider the case of large h.
As h/J increases from the value of 0.5, we observe that the form of P (s) changes from
PWig (s) to PPoi (s). The observation suggests that the eﬀect of random magnetic ﬁeld oneach site should become larger than that of the correlation between adjacent spins, as therandom ﬁeld h/J increases. The spins should become more independent of each other as
h/J increases, since the eﬀect of correlation decreases eﬀectively. Thus, the Poisson-likebehavior of P (s) should appear in the limit of large h/J. Similar shifts from PWig (s) to
PPoi (s) as randomness increases have been discussed for the 3D Anderson model [40], thespin-glass clusters [37], and the open-boundary Heisenberg chain [39].
For the case of small h, the level-spacing distribution P (s) strongly depends on thesystem size L, and the behavior of P (s) is dominated by ﬁnite-size eﬀects. In Fig. 4.11, we
observe that the form of P (s) for h/J = 0.1 is diﬀerent from that of the standard Wigner
distribution, particularly when L is small. When L is small, random magnetic ﬁeld isirrelevant to energy levels if it is smaller than the order of 1/L. In fact, energy diﬀerences
should be at least in the order of 1/L, and random magnetic ﬁeld can be neglected if it ismuch smaller than some multiple of 1/L. Thus, for the case of small h, the level statistics
should show such a behavior as that of h = 0. In fact, the Hamiltonian for h = 0 is theintegrable XXZ spin chain, which should show Poissonian behavior. Furthermore, theintegrable XXZ Hamiltonian at ∆ = 0.5 has the sl 2 loop algebra symmetry [7], and the

43

1.4

h/J=0.1
h/J=0.5
h/J=1
h/J=2

1.2

η

1

Poisson-like

0.8
0.6
0.4
0.2
0

Wigner-like
8

10

L

12

14

Figure 4.13: Dependence of the parameter η on the lattice size L for h/J =0.1, 0.5, 1,
2 and ∆ = 0.5 (reprinted form Ref. [16]). η = 0 corresponds to the Wigner distributionand η = 1 to the Poisson distribution.

level-spacing distribution should show a peak at s = 0 as shown in Fig. 4.4. In Fig. 4.11,
the graph of P (s) for L = 8 and h/J = 0.1 suggests such behavior.
The behavior of the level-spacing distributions can also be conﬁrmed in Fig. 4.12,
where ∆3 (l) are shown for the four cases corresponding to those of Fig. 4.11. The behaviorof ∆3 (l) changes from Wigner to Poissonian behavior as h/J becomes larger for eachL, though the behavior of that for h/J = 0.1 strongly depends on the system size L.
Particularly, the slope of ∆ 3(l) for L = 8 and h/J = 0.1 is larger than that of Poissonian
behavior.

4.5.3

Dependence on the system size

Let us discuss the ﬁnite-size eﬀects on the level-spacing distributions. In order to observe
the size dependence of P (s) clearly, we employ the parameter η deﬁned by Eq. (4.15). In
Fig. 4.13, the value of η for h/J = 0.1 strongly depends on the lattice size L. Moreover,
we observe that as L increases, η decreases for h/J = 0.5, while η increases for h/J = 2.
The observation suggests that η approaches either the value 0 or 1 as L increases. In
other words, it should become more deﬁnite whether P (s) has Wigner-like behavior ornot, as the system size becomes large.

44

1

(a)
0.8

P(s)

0.6
0.4
0.2
0
0

0.5

1

s

1.5

2

2.5

3

1

(b)
0.8

P(s)

0.6
0.4
0.2
0
0

0.5

1

s

1.5

2

2.5

3

Figure 4.14: Level-spacing distribution of the XXZ spin chains under random magnetic
ﬁeld for L = 14, ∆ = 0.5, and (a) h/J = 0.1; (b) h/J = 0.5. Red histograms are for all
levels; green bars, 1/3 of all levels around the center; blue crosses, 10% of all levels fromeach of the two edges. Solid and broken lines show the Wigner and Poisson distributions,
respectively.

45

4.5.4

Homogeneity of the characteristic behavior of level statistics throughout the spectrum

The characteristic behavior of level statistics does not depend on the energy range of thespectrum for XXZ chains under random magnetic ﬁeld, either. In Fig. 4.14, we show
the level-spacing distributions for (a) h/J = 0.1 and (b) h/J = 0.5 when L = 14 and∆ = 0.5. There is no remarkable diﬀerence among the diﬀerent energy ranges. Even whenh/J = 0.1, P (s) shows only little diﬀerence from PWig (s) for the energy range near theedges.

4.6

Summary

For the ﬁnite spin- 12 XXZ spin chains with the NNN interaction or with a random magnetic
ﬁeld, we have evaluated characteristic quantities of level statistics such as the level-spacing
distribution, the spectral rigidity and the number variance. We have conﬁrmed the correspondence between non-integrability and Wigner behavior in the spectrum. We havealso discussed various cases where unexpected non-Wigner behavior appears. Throughthe numerical results we have obtained the following conjecture: When the symmetryof a ﬁnite-size system enhances at some points of the parameter space, the characteristic
behavior of level statistics should be given by Poisson-like behavior near some region close
to the points. In particular, we have shown that ﬁnite-size eﬀects play an important rolein the characteristic quantities of level statistics for the XXZ spin chains.
The NNN coupled XXZ chains are integrable for J 2 = 0, and their U (1) symmetry
extends into SU (2) symmetry at the point of ∆1 = ∆2 = 1. When the parameters arenear J2 = 0 and/or ∆1 = ∆2 = 1, the behavior of level statistics is Poisson-like behavior.
We have also shown that some unexpected behavior can appear when an extra symmetry
is not considered for desymmetrization, such as the case of the spin reversal symmetry inz= 0. Here we note that in some cases extra symmetries depend on somethe sector of Stot
parameters as well as some quantum numbers.
For ∆ = 0, Poissonian behavior appears in the level statistics of the XXZ spin chains
under random magnetic ﬁeld although the system has the randomness. As ∆ increases
from zero, P (s) rapidly shifts to the Wigner distribution. The Poissonian behavior of
P (s) and ∆3(l) have been explained in terms of Anderson localization. For ∆ = 0, P (s)and ∆3 (l) strongly depends on L when h is small. When L is ﬁnite, P (s) and ∆3(l) should
show Poisson-like behavior in the small h limit. In the large h limit, however, P (s) and∆3(l) should become close to the Poissonian behavior independent of L.

46

Chapter 5
Energy Diﬀusion
In this chapter, we study the dynamical properties of energy spectra. We investigate energy diﬀusion for next-nearest-neighbor (NNN) coupled XXZ spin chains in the presenceof a periodically oscillating magnetic ﬁeld. Diﬀusion coeﬃcients are found to obey thepower law with respect to both the ﬁeld strength and driving frequency with its power
varying depending on the linear response and non-perturbative regimes [18]. The widthsof the linear response and the non-perturbative regimes depend on the strength of frustrations. We have also elucidated a mechanism for oscillation of energy diﬀusion in thecase of weakened frustrations.

5.1

Preliminaries

There exists an accumulation of studies on quantum dynamics of classically chaotic systems, e.g. kicked rotators, kicked spin-tops, hydrogen atoms in time-dependent electric
ﬁeld, and the standard map model, to mention a few [45]. Quantum suppression of energy
diﬀusion, dynamical localization and other signatures of quantum chaos are notable inthese dynamics. However, most of the systems treated so far are conﬁned to those with afew degrees-of-freedom, and little attention is paid to dynamics of quantum many-bodysystems [46, 47, 48] whose adiabatic energy levels are characterized by Gaussian orthogonal ensemble (GOE) spectral statistics, i.e. by a hallmark of quantum chaos. Whilesome important contributions [49, 50, 51, 52, 53, 54, 55] are devoted to dynamics of akind of many-body systems, those systems are actually described by the random-matrix
models, and not by deterministic quantum Hamiltonians. It is highly desirable to explore
dynamical behaviors of deterministic quantum many-body systems exhibiting GOE or
Gaussian unitary ensemble (GUE) spectral statistics.
On the other hand, the frustrated quantum spin systems have been receiving awide attention, and we can ﬁnd their realization in s = 12 antiferromagnetic chainsCu(ampy)Br2 [11] and (N2H5 )CuCl3 , [12] and in s = 12 triangular antiferromagnets [56].
47

The high-lying states of these quantum many-body systems deserve being studied in thecontext of “quantum chaos.” The advantage of the frustrated quantum systems is that
one can expect quantum chaotic behaviors appearing already in the low energy region
near the ground state [57, 58]. From the viewpoint of real physics of condensed matters, novel features observed in the low-energy region are very important and welcome.
Recalling that in most of deterministic Hamiltonian systems quantum chaotic behaviors
appear in high-lying states, the role of frustration is essential in the study of quantum
dynamics from the ground state of deterministic many-body systems with GOE or GUE
level statistics.
In this work, we investigate dynamics of XXZ quantum spin chains which have antiferromagnetic exchange interactions for the nearest-neighbor (NN) and the next-nearestneighbor (NNN) couplings. The NNN couplings cause the frustration, i.e. diﬃculty inachieving the ground state, thereby attributing a name of frustrated quantum spin chainsto these systems. In fact, the level statistics of the NNN coupled XXZ spin chains without
an applied magnetic ﬁeld has been studied intensively in Refs. [15, 17], and it has beenshown that GOE behavior, which is typical of quantum chaos, appears already in the low
energy region near the ground state.
A natural extension of the research is to investigate dynamics of the frustrated
quantum spin chains with an applied periodically oscillating magnetic ﬁeld. We calculate
a time evolution of the system starting from their ground state and analyze the natureof energy diﬀusion. We shall numerically exhibit the time dependence of energy variance,and show how the diﬀusion coeﬃcients depend on the coupling constants, the anisotropy
parameters, the magnetic ﬁeld and the frequency of the ﬁeld. Furthermore, to comparewith the energy diﬀusion in the case of weakened frustrations, we also investigate dynamics
of the corresponding energy diﬀusion in XXZ spin chains with small NNN couplings and/or
exchange anisotropy.
The organization of this chapter is as follows: In Sec. 5.2, we shortly explain theidea of energy diﬀusion. In Sec. 5.3, we brieﬂy describe a numerical approach to obtainthe time evolution operator. In Sec. 5.4, we shall show the time dependence of energy
variance starting from the ground state of the many-body system and explain a way toevaluate diﬀusion coeﬃcients. Section 5.5 elucidates how diﬀusion coeﬃcients depend on
ﬁeld strength and driving frequency. Here the universal power laws are shown to exist in
the linear response and non-perturbative regions. Section 5.6 is devoted to a mechanismof oscillation of energy diﬀusion. We summarize this chapter in Sec. 5.7.

48

E

t
Figure 5.1: Schematic picture of energy spectra.

5.2

What is energy diﬀusion like?

Let us consider a system which has a time dependent parameter X. The solution of thetime-independent Schrödinger equation gives a set of energy levels at a ﬁxed time. The
sets of energy levels arranged in time order give adiabatic energy spectra such as Fig. 5.1.
Because of the time-dependent parameter X, the system obeys the time-dependent
Schrödinger equation that causes nonadiabatic transitions among energy levels. When the
nonadiabatic transitions occur repeatedly, the occupation probabilities of energy levels
diﬀuse. The diﬀusion of occupation probabilities leads to energy diﬀusion. The behavior
of energy diﬀusion is determined by the spectral property of the adiabatic energy levels
and the changing rate Ẋ of the parameter.

5.3

Numerical procedure

We give the Hamiltonian for the NN and NNN exchange-coupled spin chain on L sites
with a time-periodic oscillating magnetic ﬁeld as
H(t) = H0 + H1 (t),

(5.1)

where
H0 = J1

Lyxz(Sjx Sj+1+ Sjy Sj+1
+ ∆Sjz Sj+1
)

j=1

+ J2

L


y
x
z
(Sjx Sj+2+ Sjy Sj+2
+ ∆Sjz Sj+2
)

j=1

−

L

j=1

Bjz (0)Sjz ,

(5.2)
49

H1 (t) =L
j=1

Bjz (0)Sjz −

L

j=1

Bjz (t)Sjz .

(5.3)

Here, S α = 12 σ α and (σ x, σ y , σ z ) are the Pauli matrices; the periodic boundary conditions(P. B. C.) are imposed. The magnetic ﬁeld B jz on jth site along the z axis is chosen toform a traveling wave:


2πj
.
(5.4)
Bjz (t) = B0 sin ωt −L
The period of Eq. (5.1) as well as Eq. (5.4) is T = 2π/ω. Because of the coexisting
spatial P. B. C., however, the eﬀective period of the adiabatic energy spectra is given by
T  = T /L = 2π/(ωL). In other words, the period of the Hamiltonian operator is T , andthe spectral ﬂow of the eigenvalues has the eﬀective period T . This periodicity property
comes from the traveling-wave form of Eq. (5.4), and is advantageous for our getting a
suﬃcient number of relevant data in each period T .
When J1 > 0 and J2 > 0, the unperturbed Hamiltonian H0 without coupling to themagnetic ﬁeld is translationally invariant and corresponds to a frustrated antiferromagnetic quantum spin model exhibiting GOE level statistics [15, 17]. If J 2 = 0 and B0 = 0, it
describes an integrable and non-frustrated model. Before calculating energy diﬀusion, we
have to consider the symmetries of the model. We divide the Hamiltonian matrix to some
sectors which have the same quantum numbers. In the Hamiltonian Eq.(5.1), total S z
z
z
(Stot
) is conserved. The eigenstates with diﬀerent S tot
are uncorrelated. For convenience,z
we choose the sector of Stot = 1. On the other hand, the non-uniform magnetic ﬁeld
breaks the translational symmetry, and leads to mixing between manifolds of diﬀerent
wave-number values.
Before proceeding to consider the time evolution of a wave function, we should note:
If we use the original Hamiltonian H(t) = H 0 + H1 (t) as it stands, the mean level spacingof eigenvalues would change depending on J2 , ∆, and B0 . To see a universal feature of theenergy diﬀusion, it is essential to scale the Hamiltonian so that the full range of adiabatic
energy eigenvalues becomes almost free from these parameters. Noting that this energy
range for the original Hamiltonian is of order of L when J1 = J2 = ∆ = 1, we deﬁne thescaled Hamiltonian H(t) = H0 + H1 (t) so that the full energy range equals L at t = 0,
which will be used throughout in the text. The Schrödinger equation is given byih̄∂|ψ(t) = H(t)|ψ(t) = [H0 + H1 (t)]|ψ(t).
∂t

(5.5)

The solution of Eq. (5.5) consists of a sequence of the inﬁnitesimal processes as|ψ(t) = U (t; t − ∆t)U (t − ∆t; t − 2∆t) · · · U (2∆t; ∆t)U (∆t; 0)|ψ(0).

(5.6)

The initial state |ψ(0) is taken to be the ground state, since our concern lies in the
dynamical behaviors starting from the many-body ground state. To calculate a time50

evolution operator U (t + ∆t; t) for each short time step ∆t, we use the fourth-orderdecomposition formula for the exponential operator [59]:U (t + ∆t; t) = S(−ip 5∆t/h̄, t5 )S(−ip4∆t/h̄, t4) · · · S(−ip1∆t/h̄, t1),
where,







(5.7)



xH1 (t)xH1 (t)S(x, t) = exp
exp(xH0 ) exp
.
2
2

(5.8)

Here, t j ’s and pj ’s are the following:tj = t + (p1 + p2 + · · · + pj−1 + pj /2)∆t,p = p1 = p2 = p4 = p5 ,= 0.4144907717943757 · · ·p3 = 1 − 4p.

(5.9)

The numerical procedure based on the above decompositions is quite eﬀective when H 1 (t)and H0 do not commute and each time step is very small. Our computation below isz= 1 manifold involves 210 levels.
concerned mainly with the system of L = 10, whose Stot
To check the validity of our assertion, some of the results will be compared to those forz= 4 whose manifold involves 364 levels.
the system of L = 14 and Stot

5.4

Time dependence of energy variance

We calculate time evolution of the state and evaluate energy variances at each integer
multiple of the eﬀective period T  = T /L = 2π/(ωL). As mentioned already, we choose
the ground state as an initial state, following the spirit of real physics of condensed
matters. This viewpoint is in contrast to that of the random matrix models where initialstates are chosen among high-lying ones [49, 50, 51, 52, 53, 54]. Consequently, the energy
variance of our primary concern is the variance around the ground state energy E0 and isdeﬁned by(5.10)
δE(t)2 = ψ(t)|[H(t) − E0 ]2|ψ(t).
Time evolution of δE(t) 2 is shown in Fig. 5.2. The parameters except for ω are ﬁxed. The
larger ω is, the faster the energy diﬀusion grows, which is consistent with our expectations.
The details are explained in Sec. 5.5. A linear growth of δE(t) 2 in time during the ﬁrst
period can be recognized more clearly for L = 14 than for L = 10, because the numberof data points per period is increased for larger L.
For wide parameter values of the next-nearest-neighbor (NNN) coupling J 2 and theanisotropy parameter ∆, the early stage of quantum dynamics becomes to show the normal
diﬀusion in energy space, i.e. a linear growth of δE(t) 2 in time. While we investigate this
51

18
16
14

(a) L=10
ω=0.5
ω=0.4
ω=0.3
ω=0.2
ω=0.1

30
25
20
2

10

δE

2

12
δE

(b) L=14

35

8

15

6

10

4
5

2
0

0 0.5 1 1.5 2 2.5 3time (t/T)

0

0 0.5 1 1.5 2 2.5 3time (t/T)

Figure 5.2: Time evolution of energy diﬀusion for (a) L = 10 and (b) L = 14 (reprinted
from Ref. [18]). The unit of time is the period T . The parameters are the following:
J1 = J2 = 1.0, ∆ = 0.3, B0 = 1.0.

normal diﬀusion process, energy variances will ﬁnally saturate because the system size
we consider is ﬁnite. On the other hand, energy variances can also saturate because ofanother reason, i.e. the dynamical localization eﬀect. It is associated with a periodic
perturbation. In any case, diﬀusion coeﬃcients have to be determined for times where
saturation does not yet occur. We determine the diﬀusion coeﬃcient D from the ﬁtting
δE(t)2 = Dt + const.

(5.11)

to some data points around the largest slope in the ﬁrst period, where the normal diﬀusion
is expected.

5.5

Diﬀusion coeﬃcients: dependence on ﬁeld strengthand frequency

Since the time evolution of our system starts from the ground state, we consider nonadiabatic regions where inter-level transitions frequently occur. In other words, we suppress a near-adiabatic or the so-called Landau-Zener (LZ) region where the driving frequency ω is much smaller than the mean level spacing divided by Planck constant. Becauseof a large energy gap between the ground and ﬁrst excited states, the near-adiabatic region cannot result in the notable energy diﬀusion and will be left outside a scope of thepresent study.

52

(a) J2=1.0(b) J2=0.2

100

100

10

10
D

D
1

1

0.1

0.01

0.5

1

ω

B0=2.0B0=1.5B0=1.0B0=0.8B0=0.55B0=2.0B0=1.5B0=1.0B0=0.8B0=0.5

0.1

0.01

0.5

1

ω

5

Figure 5.3: Driving frequency dependence of the diﬀusion coeﬃcients (reprinted from
Ref. [18]). The chained line and the solid line are just eye guides for D ∝ ω β with
β = 1 and 2, respectively. The symbols () are the average of the diﬀusion coeﬃcients
calculated for several values of ∆ (0.3 ≤ ∆ ≤ 0.8). The parameters are the following:
L = 10, J1 = 1.0; (a) J2 = 1.0, (b) J2 = 0.2.

Beyond the LZ region, however, so long as the changing rate Ẋ of a perturbation
parameter is not very large,1 the diﬀusion coeﬃcient can be calculated using the Kubo
formula. We call such a parameter regime “linear response” regime. In the linear response regime, D ∝ Ẋ 2 (See, e.g., Refs. [50] and [51]). When Ẋ is large, however, theperturbation theory fails. We call such a parameter regime “non-perturbative” regime.
In the non-perturbative regime, the diﬀusion coeﬃcient is smaller than that predictedby the Kubo formula [51, 53]. According to Ref. [51], D ∝ Ẋ γ with γ ≤ 1 in the nonperturbative regime. We note that Ẋ ∝ B0 ω in this work since the perturbation is given
by Eq. (5.4). Both Refs. [51] and [53] are based on the random matrix models, which areutterly diﬀerent from our deterministic one. The discussion about energy diﬀusion for arandom matrix model [51] is brieﬂy reviewed in Appendix E.
Numerical results of diﬀusion coeﬃcients in Fig. 5.3 are almost consistent with theargument of Ref. [51]. Diﬀusion coeﬃcients as a function of ω are shown in Fig. 5.3. In
Fig. 5.3(a), where J2 = 1.0 (i.e. the fully-frustrated case), D is larger as B 0 is larger fora ﬁxed value of ω. In a small-ω regime, D ∝ ω β with β = 2, though β > 2 for small B0 .
The latter is merely attributed to the fact that the perturbation is too small to observea suﬃcient energy diﬀusion when both ω and B0 are small. In a large-ω regime, β = 1.
1

More precisely, this phrase means that the time scale upon which the time-dependent Hamiltonian
matrix elements decorrelate is much larger than the time scale corresponding to the typical separation ofenergy levels.

53

60
ω=3
ω=5

40

δE

2

50

30
20
10
0

0

1

2
time (t/T)

3

Figure 5.4: Time evolution of energy variance for large ω. The parameters are L = 10,J1 = J2 = 1.0, ∆ = 0.3, B0 = 2.0.

Namely, we observe that β = 2 in the linear response regime and β = 1 in the nonperturbative regime. In fact, for a large-ω regime, the increase of energy variances per
eﬀective period hardly depend on ω until δE(t) 2 starts to decrease as shown in Fig. 5.4.
This explains the observation that D ∝ ωβ with β = 1 in both Fig. 5.3(a) and Fig. 5.3(b).
Let us represent the increase of energy variances per eﬀective period as ∆(δE 2). From
the deﬁnition of D, i.e. Eq. (5.11), D ∝ ∆(δE 2)/T . Then we ﬁnd that D ∝ ω if ∆(δE 2)
is constant.
On the other hand, in Fig. 5.3(b) where J2 = 0.2 (i.e. a weakly-frustrated case), theregion with β = 1 expands. For small B0 , β > 2 in a small-ω regime. For small B 0 andaround ω ∼ 1, D seems to decrease rather than increase especially in the case of J 2 = 0.2.
In fact, the slope of energy variance for large ω is larger than that for small ω in the very
early stage of energy diﬀusion. However, some kind of localization occurs for large ω asshown in Fig. 5.5, which leads to the suppression of D.
It is seen more clearly in Fig. 5.6 how the behavior of D changes between a linear
response regime and a non-perturbative regime. The diﬀusion coeﬃcient D obeys thepower law D ∝ (B0ω)β with its power β being two in the linear response regime andβ = 1 in the non-perturbative regime. For small B 0 ω, the power law seems to fail becauseof some ﬁnite-size eﬀects. These universal feature is conﬁrmed in systems of larger size.
Actually, D obeys the power law better for L = 14 [Fig. 5.6(b)] than L = 10 [Fig. 5.6(a)].
In addition, error bars are shorter for L = 14 than L = 10. Here, we have used the dataof ω ≤ 1. We cannot expect meaningful results in a large-ω regime since, as mentioned
above, energy diﬀusion is not normal there.
Figure 5.6 suggests that the strength of frustration should aﬀect the range of the
54

35
30

ω=1
ω=1.5
ω=2

10

0.4

8

25

6

20

4

15

2
0

2

δE

0.2

10
5
0
0

0.5

1

1.5
time (t/T)

2

2.5

3

Figure 5.5: Time evolution of energy variance for small B 0 and around ω ∼ 1. The
parameters are L = 10, J1 = 1.0, J2 = 0.2, ∆ = 0.3, B0 = 0.8. The inset is the magniﬁed
ﬁgure of the hatched part.

linear response regime. The linear response regime is shorter for J 2 = 0.2 than for J2 = 1.0,
while the non-perturbative regime is larger for J 2 = 0.2 than for J2 = 1.0. In fact, whenJ2 = 0 (i.e. the integrable case), D ∝ (B 0ω)β with β = 1 for almost all the data in thesame range of B0ω as that of Fig. 5.6.

5.6

Oscillation of energy diﬀusion in weakly-frustratedcases

We shall now proceed to investigate oscillations of diﬀusion which occur in the nonperturbative regime of a weakly-frustrated case. Figure 5.7(a) shows an example of oscillatory diﬀusion for J2 = 0.2, which is compared with a non-oscillatory diﬀusion for
J2 = 1.0. The two examples have the same set of parameters except for J 2. However,the cases of J2 = 1.0 and J2 = 0.2 are in the linear response regime and in the nonperturbative regime, respectively. The variance for both cases shows normal diﬀusion at
the very early stage of time evolution. For J 2 = 1.0, the energy variance seems to saturate
after a normal diﬀusion time. On the contrary, the energy variance for J 2 = 0.2 shows
large-amplitude oscillations. To investigate more details, we introduce another deﬁnitionof energy variance:
δ Ẽ(t)2 = ψ(t)|[H(t) − ψ(t)|H(t)|ψ(t)] 2|ψ(t).

55(5.12)

100

100(b) L=14(a) L=10
10

10

1

1

D

D

0.1

0.1

10

10
1

1

0.01

0.01

0.1

0.1
0.3

0.1

B0ω

1

1

0.3

3

5

0.1

B0ω

1

1

3

5

Figure 5.6: Dependence of the diﬀusion coeﬃcients on the product of ﬁeld strength B 0
and driving frequency ω for (a) L = 10 and (b) L = 14 (reprinted from Ref. [18]). The
symbols () are the average of the diﬀusion coeﬃcient calculated for several values of ∆(0.3 ≤ ∆ ≤ 0.8). For the main panels, J 1 = J2 = 1.0; for the inset, J1 = 1.0 and J2 = 0.2.
The chained line and the solid line are just eye guides for D ∝ (B 0 ω)β with β = 1 and 2,
respectively. Some error bars are too short to see.

56

40

(a)

35

δE

2

30
25
20
15
10
J2=1.0J2=0.2

5
0

0

0.5

1

1.5
time (t/T)

2

2.5

3

0.5

1

1.5
time (t/T)

2

2.5

3

6

~2
δE

(b)
4
2
0

0

Figure 5.7: Examples for time evolution of energy variances : (a) δE(t) 2 and (b) δ Ẽ(t)2(see text) (reprinted from Ref. [18]). Solid lines are for J 2 = 1.0; Broken lines, J 2 = 0.2.
The parameters are the following: L = 10, J1 = 1.0, ∆ = 0.3, B0 = 1.5, ω = 0.5.

57This follows a standard deﬁnition of the variance and quantiﬁes the degree of energy
diﬀusion around the time-dependent expectation of the energy Hamiltonian. The time
evolutions of δ Ẽ(t)2 corresponding to that of δE(t)2 are shown in Fig. 5.7(b). In the
fully-frustrated case (J2 = 1.0), the proﬁle of δ Ẽ(t)2 is similar to that of δE(t)2. Thisobservation indicates that an occupation probability spread over the whole levels after
normal diﬀusion of energy.
On the contrary, in a weakly-frustrated case (J2 = 0.2) in Fig. 5.7, δ Ẽ(t)2 shows
small-amplitude oscillations reﬂecting the large-amplitude oscillations of δE(t) 2. Most
part of δ Ẽ(t)2 for J2 = 0.2 is smaller than that for J2 = 1.0. Furthermore, minima ofδ Ẽ(t)2 come just before minima and maxima of δE(t) 2. These observations indicatesthe following: an occupation probability, which is diﬀusing slowly, clustering aroundthe expectation of energy oscillates together with the expectation in the energy space.
To make the picture of such behavior clearer, let us consider an occupation probability
described by(5.13)Pt (En ) = | φn |ψ(t)|2,
where |φn  is the nth excited eigenstate of H 0 :H0 |φn  = En |φn .

(5.14)

When t = 0 , Pt (En ) is given by the Kronecker delta: P 0 (En ) = δEn,E0 , where E0 is theenergy of the ground state at t = 0. As t increases, Pt (En ) forms a wave packet in energy
space and moves to higher levels. When the wave packet reaches some highest levels, it
reﬂects like a soliton and moves back to lower levels. Such behavior is repeated, although
the wave packet of Pt (En ) broadens slowly. We have actually watched this soliton-likebehavior of Pt (En ) in a form of an animation.
The picture discussed above is also supported by the adiabatic energy spectra in
Fig. 5.8. Figures 5.8(a) and 5.8(b) correspond to fully- and weakly-frustrated cases,
respectively. Much more sharp avoided crossings appear in Fig. 5.8(b) than Fig. 5.8(a).
Some energy levels appear to be crossing, although they are very close and never crossing
in fact. At a sharp-avoided-crossing point, Landau-Zener formula for two adjacent levels
is applicable. Then the nonadiabatic transition leads to one-way transfer of a populationfrom a level to its partner, failing to result in the energy diﬀusion. For small-J 2, therefore,
successive sharp avoided crossings can suppress diﬀusion of energy.
We believe that large-amplitude oscillations of δE(t) 2 should be one of characteristicfeatures of the non-perturbative regime in this ﬁnite frustrated spin system. In fact,
similar oscillations of energy variance are seen for large ω and large B0 even when J2 = 1.0
though the energy variance rapidly converges after one or two periods. For example, sucha behavior of energy variance is seen in Fig. 5.4. How long such oscillations continue
should depend mainly on J2 .
58

(a) J2=1.0

1.2(b) J2=0.2

1.2

E

E
1

0.8

1

0

π/4

ωt

π/2

0.8

π/4

0

ωt

π/2

Figure 5.8: Parts of energy spectra depending on adiabatically ﬁxed time t with 0 ≤ t ≤T /4 (reprinted from Ref. [18]). Eﬀective period is ωT  = 2π/10. The parameters are thefollowing: L = 10, J1 = 1.0, ∆ = 0.3, B0 = 0.8; (a) J2 = 1.0, (b) J2 = 0.2.

1
1

2

0.8

0.8
0.6

P(s)

0.4

0.6

0.2

0.4
0.2
0
0

0.5

1

s

1.5

2

2.5

3

Figure 5.9: Level-spacing distributions at t = π/4 for lowest 300 levels from the groundstate (about 10% of all 3003 levels) (reprinted from Ref. [18]). Blue histogram is for
J2 = 1.0; Red bars, J2 = 0.2; Solid curve, GOE spectral statistics. The other parametersz= 1, J1 = 1.0, ∆ = 0.3, B0 = 0.8. The inset is for all
are the following: L = 14, Stotlevels when J 2 = 1.0. The numerical methods to obtain the level-spacing distributions
are referred in Refs. [15, 17].

59

It is a notable fact that, common to both J2 = 1.0 and J2 = 0.2, the level-spacing
distributions in Fig. 5.9 show GOE behavior. This GOE behavior in the adiabatic energy
spectra appears for an arbitrary ﬁxed time except for special points such as t = T = 2π/ω.
This fact suggests that dynamics can reveal some various generic features of quantum
many-body systems which can never be explained by level statistics. The level-spacing
distributions in Fig. 5.9 convey another crucial fact: they have been calculated for low
energy levels because our interest is in the low energy region around the ground state.
We have conﬁrmed that the level-spacing distributions for all energy levels in the inset
is also described by GOE spectral statistics. It is typical of this frustrated spin systemthat GOE level statistics is observed already in the low energy region [17]. We havealready discussed homogeneity of the behavior of level statistics throughout the spectrumin Chapter 4. One may wonder why P (s) for J2 = 0.2 shows Poisson-like behavior in
Fig. 4.9 (b) and GOE behavior in Fig. 5.9. This is because a magnetic ﬁeld is added inthis model, Eq. (5.1).

5.7

Summary

We have explored the energy diﬀusion from the ground state in frustrated quantum XXZ
spin chains under the applied oscillating magnetic ﬁeld. In a wide parameter region of
next-nearest-neighbor (NNN) coupling J2 and exchange anisotropy ∆, the diﬀusion isnormal until the dynamical saturation begins to prevail. Diﬀusion coeﬃcients D obeys
the power law with respect to both the ﬁeld strength and driving frequency with its
power being two in the linear response regime and equal to unity in the non-perturbative
regime. In the case of weakened frustrations with small-J 2 we ﬁnd oscillation of energy
diﬀusion, which is attributed to a non-diﬀusive and ballistic nature of the underlying
energy diﬀusion. In this way, the energy diﬀusion reveals generic features of the frustrated
quantum spin chains, which cannot be captured by the analysis of level statistics.

60

Chapter 6
Conclusions
In this thesis, we have discussed the properties in the energy spectra of XXZ spin chains.
The XXZ spin chains are basic quantum spin chains and related to the most importantfamilies of integrable spin systems. We have studied level crossings, level bifurcations,and level statistics for the statistical properties, and energy diﬀusion for the dynamical
properties.
To study level crossings, the eigenvalues has been calculated for 3-site, 4-site and6-site XXZ spin chains. We have observed many counterexamples to the noncrossing rulein the spectra of the integrable XXZ spin chains. Several level crossings are related tothe sl2 loop algebra symmetry. The symmetry of a XXZ spin chain enhances at someparticular values of the anisotropy parameter, where multiple level crossings are observed.
Branches have been observed in the ﬂow of eigenvalues for the inhomogeneous transfer matrix of the XXZ spin chain. Those branches can be called “level bifurcations”. The
appearance of the branches has been conﬁrmed through the Bethe ansatz method. The
observation of level bifurcations may contribute to the development on the exact calculation of thermodynamic quantities for the XXZ spin chain. The transfer matrix of a spin
chain is often used to calculate the thermodynamical quantities of it.
Level statistics has been studied for XXZ spin chains with next-nearest-neighbor(NNN) coupling or with a random magnetic ﬁeld. We have conﬁrmed correspondence
between non-integrability and Wigner behavior in the energy spectra and shown how
ﬁnite-size eﬀects aﬀect the level statistics of the ﬁnite-size XXZ spin chains. Studying theﬁnite-size eﬀects on the level statistics should be useful for investigating similar systems
or for calculating quantities such as energy diﬀusion. It has also been suggested that
incomplete desymmetrization can cause unexpected behavior in level statistics. Study ofthe unexpected behavior should be helpful to resolve similar behavior of level statistics
even in other systems.
Energy diﬀusion has been investigated to study the dynamical properties for NNN
coupled XXZ spin chains under a periodically oscillating magnetic ﬁeld. We have observed61that diﬀusion coeﬃcients D shows a universal behavior: D obeys the power law withrespect to both the ﬁeld strength and driving frequency. The observation may lead tosome new nature of frustrated quantum spin systems. For weakened frustrations, wehave found the existence of small ﬂuctuation of diﬀusion coeﬃcients and the non-diﬀusive
behavior of the energy diﬀusion. The behavior should suggest that energy diﬀusion reveals
some features, which cannot be captured by the analysis of level statistics for the XXZ
spin chains.
In this thesis, we have focused on the properties in the energy spectra of XXZ spin
chains. We believe that the above intensive study on the XXZ spin chains is very useful forinvestigating other ﬁnite quantum spin chains or ﬁnite low-dimensional quantum systems.

62

Acknowledgements
I wish to thank my supervisor, Professor Tetsuo Deguchi, for his advice, guidance, discussion, concern and encouragement.
I would like to thank Professor Fumiaki Shibata and Professor Katsuyoshi Kobayashi
for useful comments and suggestion.
I am grateful to Professor Katsuhiro Nakamura of Osaka City University for collaboration in the work of Chapter 5 and for useful discussion on the work of Chapter 4. I also
wish to acknowledge valuable discussion on the work of Chapter 4 with Professor Takeo
Kato of University of Tokyo. I also thank to the stuﬀs and students of Prof. Nakamura’sLaboratory at Osaka City University.
I would like to express my special thanks to my friends whom I met at Bussei-Wakate
Summer School and Butsurikko-Club Seminar. In particular, I would like to thank Mr.
Hirofumi Watanabe for telling me how to perform unfolding.
I wish to thank Professor Yshai Avishai of Ben-Gurion University in Israel for valuable comments and suggestion.
I wish to acknowledge the Yukawa Institute for Theoretical Physics at Kyoto University. Discussions during the YITP workshop YITP-W-03-16 on “Quantum Mechanics
and Chaos: From Fundamental Problems through Nanosciences” were useful to complete
the work of Chapter 4.
This study is partially supported by Hayashi Memorial Foundation for Female Natural Scientists.
Finally, I wish to thank the members of Deguchi Laboratory, my friends and myfamily for support and encouragement.

63

Appendix A
Hamiltonian Matrix
Let us show how to describe Eq. (2.1) in the form of the Hamiltonian matrix. The matrices
on the jth site, σj± and σjz , are given by· · ⊗ I$ ⊗ !"#$σ α ⊗ I! ⊗ ·"#· · ⊗ I$ .
σjα = I! ⊗ ·"#
jth

j−1

(A.1)

N −j

Here,




0 1 
,
σ+ = 
0 0





0 0 
σ− = 
,
1 0





1 0 
σz = 
,
0 −1





1 0 
I=
.
0 1

(A.2)

Now we remember








a11 a12   b11
A⊗B = ⊗a21 a22b21

a11 b11b12   a11 b21
=
 a b
b22
 21 11
a21 b21a11 b12a11 b22a21 b12a21 b22

a12 b11a12 b21a22 b11a22 b21

a12 b12
a12 b22 
a22 b12 a22 b22(A.3)

For example, when L = 3, Eq. (2.1) is described by




3∆
0
0
0
0
0
0
0



0
0 
−∆
2
0
2eiφ 0
 0



 0
2
−∆
0
2
0
0
0 




iφ

0 
0
0
−∆
0
2 2e
J 0

H= 
.
−iφ
4
2
0
−∆ 0
0
0 
 0 2e



 0
0
0
2
0 −∆ 2
0 




−iφ
 0

0
2
−∆
0
0
0
2e


0
0
0
0
0
0
0 3∆

64

(A.4)

Here, we give the base vectors by | ↑↑↑, | ↑↑↓, | ↑↓↑, · · ·, | ↓↓↓. Namely,












1




 0 




 0 
1



 
 





 



1   1   1   0   1   0 | ↑↑↑ =⊗⊗= ⊗=  ,
 0 
 0 
0
0
0
0








0
 0 


 0 


0
0


 1 








1
0



 
 





 




1
1
0
0
0
0




| ↑↑↓ =   ⊗   ⊗   =   ⊗   =   ,




0
0
1
1
 0 
 0 


 0 
0




 0 


0
0




 0 




 1 
0



 
 





 




1
0
1
1
1
0




| ↑↓↑ =   ⊗   ⊗   =   ⊗   =   , · · · .




0
1
0
0
 0 
 0 


 0 
0




 0 


0

65

(A.5)

Appendix B
Jordan-Wigner and Fourier
Transformations
Let us rewrite Eq. (4.3) by
H = H1 + H2
= J1
+ J2

L

l=1
L

l=1

y
x
z
(Slx Sl+1+ Sly Sl+1
+ ∆1Slz Sl+1)yxz(Slx Sl+2+ Sly Sl+2
+ ∆2Slz Sl+2
),

(B.1)

where H1 is the term containing nearest-neighbor couplings and H2 is the term containing
next-nearest-neighbor couplings. We deﬁne the Jordan-Wigner transformation by

σl− = exp −iπ


l−1

†
j=1



cj cj  c†l ,

(B.2)


l−1

exp iπc†j cj  cl ,(B.3)

σlz = 2σl+ σl− − 1 = 1 − 2σl− σl+ ,

(B.4)

σl+ =

j=1

where c†j and cj are the creation and annihilation operators of fermions on jth site. Under
the Jordan-Wigner transformation, H1 is written by


H1

L

J1
∆1 L − 2=c†l cl4
l=1J1 L−1[2(c†l cl+1 + c†l+1 cl ) + ∆1(4c†l cl c†l+1 cl+1 − 2c†l+1 cl+1 )]
4 l=1
J1
+
[2(c†L cL+1 + c†L+1 cL ) + ∆1(4c†L cL c†L+1 cL+1 − 2c†L+1 cL+1 )].
4

+

66

(B.5)

±
Here, we consider periodic boundary conditions (σ L+1= σ1± ):









L
iπL iπ += exp −iπ
σj− σj+  σL+1
= exp −σjz  σ1+
+
2
2
j=1
j=1

cL+1

L


'

=

σ1+ exp




(

L
z
iπ − + Stot
+1
2 

,





L
iπL iπ −= exp iπ
σj− σj+  σL+1
= exp σjz  σ1−
−
2
2
j=1
j=1

c†L+1

L


'

= σ1− exp iπ



(

L
z
+1
− Stot
2

,

(B.6)

z
= (L − M )/2 − M/2 =where we use eS1 σ1± = σ1± eS1 ±1 . Since σ 1+ = c1, σ1− = c†1 , and Stot
L/2 − M , where M is the number of fermions, we havezzc†L+1 = −(−1)M c†1 .

cL+1 = −(−1)M c1 ,Therefore, H1 is rewritten by
H1

L

J1
∆1 L + (4c†l cl c†l+1 cl+1 − 4c†l cl )
=
4
l=1

(B.7)




L−1
 †
J1
+
(cl cl+1 + c†l+1 cl ) − (−1)M (c†L c1 + c†1cL ) .
·2
4
l=1

(B.8)

Now, we deﬁne the Fourier transformation as
1  iklĉk e ,
cl = √L k
1  † −iklĉk e .
c†l = √
L k

(B.9)
(B.10)

Here, k takes (2π/L) × (an integer) for odd M and (2π/L) × (a half-integer) for even M ,and 0 ≤ k < 2π. After the Fourier transformation,
J1∆1 L + J1 (cos k − ∆1 )ĉ†k ĉk4kJ1 ∆ 1 −δk +k ,k +k exp [−i(k2 − k4 )] ĉ†k1 ĉ†k2 ĉk3 ĉk4 ,L k1 ,k2 ,k3 ,k4 1 2 3 4

H1 =

where
δk1 +k2 ,k3 +k4 =


 1
 0

when k1 + k2 (mod 2π)=k3 + k4 (mod 2π).
otherwise.

(B.11)(B.12)

Considering combination of k’s, we can rewrite Eq. (B.11) by
H1 = J1 ∆1
−






L
cos k · ĉ†k ĉk− M + J1
4k

2J1 ∆1 δk +k ,k +k [cos (k2 − k4 ) − cos (k2 − k3 )] ĉ†k1 ĉ†k2 ĉk3 ĉk4 . (B.13)
L k1 <k2 1 2 3 4
k3 >k467

In the same way, H2 can be rewritten by
H2 = J2 ∆2
+

−






L
cos 2k · ĉ†k ĉk− M + J24k

2J2 δk +k ,k +k [cos(k1 + k3 ) + cos(k2 + k4 )L k1 <k2 1 2 3 4
k3 >k4− cos(k1 + k4 ) − cos(k2 + k3 )] ĉ†k1 ĉ†k2 ĉk3 ĉk4

2J2 ∆2 δk1 +k2 ,k3 +k4 {cos [2 (k2 − k4 )] − cos [2 (k2 − k3 )]} ĉ†k1 ĉ†k2 ĉk3 ĉk4 . (B.14)
L k1 <k2
k3 >k468

Appendix C
Spin reversal symmetry on
momentum-based fermions
Let us ﬁnd the momentum-based expression of the mapping corresponding to the spinreversal transformation (Sj± → Sj∓ , Sjz → −Sjz ). According to Eqs. (B.2), (B.3), (B.9),
and (B.10),




L
l−1

1 
exp(−ikl) exp −iπ
σj− σj+  σl+ ,ĉk = √
L l=1
j=1

1
ĉ†k = √

L


L l=1



exp(ikl) exp iπ

l−1

j=1

(C.1)



σj− σj+  σl− .

(C.2)

Under the transformation σl± → σl∓ , Eq. (C.1) is transformed as the following.


ĉk



L
l−1

1 
→ √
exp(−ikl) exp −iπ
σj+ σj−  σl−
L l=1
j=1





=

L
l−1

1 
√
exp(−ikl) exp −iπ (Ij − σj− σj+ ) σl−
L l=1
j=1

=

L
l−1

1 
√
exp[−ikl − iπ(l − 1)] exp iπ
σj− σj+  σl− ,
L l=1
j=1





(C.3)

where Ij is the unit matrix. Now, considering e 2πil = 1, where l is an integer, we can ﬁndthat
exp[−ikl − iπ(l − 1)] = exp[−ikl − iπ(l − 1) + 2πil] = exp[i(π − k)l + iπ]= − exp[i(π − k)l].

(C.4)

Therefore, Eq. (C.3) is rewritten as




L
l−1

1 
exp[i(π − k)l] exp iπ
σj− σj+  σl− = −ĉ†π−k .
ĉk → − √
L l=1
j=1

69

(C.5)

In the same way, we can showĉ†k → −ĉπ−k .

70

(C.6)

Appendix D
Spin reversal operation on thevacuum state
Let us introduce the following matrix:
Us =

L

j=1

σjx .

(D.1)

We may express the spin reversal operation Eq. (4.5) as follows
Us Sj± Us−1 = Sj∓ ,Us Sjz Us−1 = −Sjz

Us ĉ†k Us−1 = −ĉπ−k ,We thus have
Us |0 =L
j=1

Us ĉk Us−1 = −ĉ†π−k .

(D.2)

σjx |0 = σ1− σ2− · · · σL− |0(D.3)

Applying the Jordan-Wigner transformation and substituting c† ’s with ĉ†k ’s through Eqs. (B.9)
and (B.10), we haveσ1− σ2− · · · σL− |0 = c†1 c†2 · · · c†L |0

1 
=
·
·
·
e−i(1k1 +2k2 +···+LkL ) ĉ†k1 ĉ†k2 · · · ĉ†kL |0
LL/2 k1kL=

1   (−i )Lj=1 jkP j )  † †e
P ĉq ĉq · · · ĉ†q |0

LL/2

1

P ∈SL

2

L

(D.4)

Thus we have the expression Eq. (4.8) of the phase factor AL. Here we recall thatqj = (2π/L)j for j = 1, 2, . . . , L for odd M , and qj = (2π/L)(j − 1/2) for j = 1, 2, . . . , L
for even M . We also recall that M denotes the number of down-spins in the sector.
We now calculate the expression Eq. (4.9) for the phase factor AL . We take a vectorz= 0 as follows. When M is odd, we introduce  = (M − 1)/2 and we deﬁne|v with Stot|v by |v = ĉ†1 ĉ†2 · · · ĉ† · ĉ†− · · · ĉ†−2 ĉ†−1 · ĉ†0 |0
71Here ĉ†j denotes ĉ†k with k = (2π/L)j. When M is even, we take  = M/2, and we deﬁne|v by |v = ĉ†1/2 ĉ†3/2 · · · ĉ†−1/2 · ĉ†−(−1/2) · · · ĉ†−3/2 ĉ†−1/2 |0
Here ĉ†j+1/2 denotes ĉ†k with k = (2π/L)(j + 1/2). Through the operation Eq. (4.6) we canshow that Us |v = (−1)+1 AL|v for M odd and Us |v = (−1) AL|v for M even. Thus,
we have at least Us |v = ±|v.
Let us show that Us |v = +|v for both odd M and even M cases. First we notethat Us c†1c†2 · · · c†L/2 |0 = c†L/2+1 c†L/2+2 · · · c†L |0. Second, expanding the vector |v in termsof c†j1 c†j2 · · · c†jL/2 |0 with 1 ≤ j1 < j2 < · · · < jL/2 ≤ L, we can show that the coeﬃcientof c†1 c†2 · · · c†L/2 |0 in the expansion is equal to that of c†L/2+1 c†L/2+2 · · · c†L |0. Therefore wehave Us |v = +|v.
Thus, we obtain the expression Eq. (4.9) for the phase factor AL.

72

Appendix E
Energy diﬀusion for a randommatrix model
Here, we review the discussion of energy diﬀusion according to Ref. [51]. First of all, we
expand the wave function in an adiabatic basis:|ψ(t) =nan (t) exp[−iθ n(t)]|φn (t),

(E.1)



1 t 
dt En (t).
(E.2)
θn (t) =h̄ 0Here, |φ n (t) and En (t) are the eigenstates and the eigenvalues of the instantaneous Hamiltonian, respectively:(E.3)
H(t)|φn (t) = En (t)|φn (t).
Substituting Eq. (E.1) into the time-dependent Schrödinger equation and using Eq. (E.2),
we obtain the equation for an (t),ȧn = Ẋφn | ∂H|φm ∂X
exp[i(θ n − θm )]am .
m=n Em − En


(E.4)

To analyze Eq. (E.4), we specify some information about the energies E n (t) and thematrix elements (∂H/∂X) nm (t) = φn |(∂H/∂X)|φm . The oﬀ-diagonal matrix elements
are Gaussian distributed with mean value zero and with a variance σ 2 which is deﬁned by

1    ∂Hσ 2(E, ∆E) = 2
ρ (E) n m=n  ∂X


nm

2

1

 δ [E − (En + Em )]δ [∆E − (En − Em )], (E.5)
2where the density of states isρ(E) =

nδ (E − En ).

73

(E.6)

Here, E = 12 (En + Em ) and ∆E = En − Em ; δ (x) is a smoothed delta function, e.g. the√Gaussian function: exp(−x2 /22 )/( 2π). We characterize the energy scale in σ 2(E, ∆E)
by the parameter ∆E0 which is deﬁned by
1∆E0 = 2σ0

 ∞0d∆E σ 2(EF , ∆E),

(E.7)

where σ02 = σ 2(EF , 0). The function σ 2(E, ∆E) in the variable ∆E decays rapidly as∆E → ±∞. The parameter ∆E 0 indicates the size of the interval over which σ 2 issigniﬁcantly diﬀerent from zero.
We assume that the matrix elements have a correlation function C(t):


∂H
∂X





∂H
(t)
∂X
nm

∗




n m

(t ) = δnn δmm σ 2(E, ∆E) C(t − t ),(E.8)

where n ≥ m and n ≥ m. The characteristic time scale for the decay of correlations of
the matrix elements is τ c = (σ0 ρẊ)−1 .
Here we consider a further simpliﬁed model since Eq. (E.4) is hard to solve directly.
We replace Eq. (E.4) byȧn (t) = Znm (t) exp[iλ(n − m)t] a m(t),

(E.9)

m=n
∗
, Znn = 0. The real
where Znm (t) are elements of an anti-Hermitian matrix: Z nm = −Zmn
and imaginary parts of Znm (t) are independently Gaussian distributed random numbers.
The correlation function of Znm (t) is

Znm (t1)Zn∗ m (t2) = δnn δmm f(n − m) C(t1 − t2 ),(E.10)

where C(0) = 1. This is a reasonable model for Eq. (E.4) when  = Ẋ , λ = (h̄ρ)−1 , andf(∆n) = σ 2(∆E)/∆E 2.
The solution of Eq. (E.9) isan (t) =m

Unm (t, 0) am (0),(E.11)

where the matrix elements of the evolution operator Û = {Unm } satisfy a Dyson equation:
Unm (t, 0) = δnm + 

 t
k

0

dt Znk (t)Ukm (t, 0) exp[iλ(n − k)t].

(E.12)

Now we derive the equation of motion for occupation probabilities: P n (t) = |an (t)|2.
Substituting Eq. (E.12) into Eq. (E.11) and using some approximation, we can obtain thefollowing (See Ref. [51] for details):
Pn (t) − Pn (0) ∼  t
2

 ∞
m

−∞

∗
dτ Znm (τ )Znm(0) Umm (τ ) exp[iλ(n − m)τ ](Pm (0) − Pn (0)).

(E.13)
74

From Eq. (E.13),

dPn 
Rnm (Pm − Pn ).
=dtm(E.14)

Here the rate constants Rnm are given by
Rnm = 2

 ∞

−∞

∗
dτ exp[iλ(n − m)τ ] Znm(τ )Znm(0) Umm (τ )

= 2f (n − m) ∞dτ exp[iλ(n − m)τ ]C(τ )c(τ ),−∞(E.15)

where c(τ ) = Umm (τ ). If the averaged occupation probability varies slowly as a functionof n, Pm − Pn can be Taylor expanded in ∆n = m − n:Pm = Pn +∂P
1 ∂ 2P∆n2 + O(∆n3 ).
∆n +
∂n
2 ∂n2

(E.16)

Therefore, Eq. (E.13) can be approximated by a diﬀusion equation:
∂ 2P
∂P
=D 2,
∂t
∂n

(E.17)

1
Rnm (n − m)2.
2 m

(E.18)

where the diﬀusion coeﬃcient isD=

To analyze Eq. (E.4), we rewrite Eq. (E.15) as
Rnm = Ẋ 2

σ(∆E)
∆E 2

 ∞
−∞

dτ exp[i∆Eτ /h̄]C(τ )c(τ ),(E.19)

where ∆E = En − Em  (n − m)/ρ. From Eq. (E.18), the diﬀusion coeﬃcient is

∞
1 3 2 ∞ρ Ẋd∆E σ 2(∆E)
dτ exp[i∆Eτ /h̄]C(τ )c(τ )
D =
2
−∞
−∞
 ∞
1 3 2
ρ Ẋ
d∆E σ 2(∆E)F (∆E),=
2−∞(E.20)

where we approximate the summation in Eq. (E.18) by an integral. The total weight of
the function F (∆E) is ∞−∞d∆E F (∆E) =

 ∞

−∞

= 2πh̄dτ C(τ )c(τ )
 ∞
−∞

 ∞

−∞

d∆E exp[i∆Eτ /h̄]dτ δ(τ )C(τ )c(τ ) = 2πh̄C(0)c(0) = 2πh̄.

(E.21)

h̄ρσ0Ẋ, F (∆E) decays rapidly to zero while σ(∆E) is signiﬁcantly diﬀerent
When ∆E0
from zero. Then the diﬀusion coeﬃcient is approximated as follows:

1 3 2 ∞ρ Ẋd∆E σ 2(∆E) F (∆E)
D =
2
−∞

1 3 2 2 ∞
ρ Ẋ σ0d∆E F (∆E)∼
2−∞
= πh̄ρ3 σ02Ẋ 2 .
75(E.22)

When ∆E0  h̄ρσ0Ẋ, σ(∆E) decays rapidly to zero while F (∆E) is signiﬁcantly diﬀerent
from zero. Then,


∞
1 3 2
ρ Ẋ τcd∆E σ 2 (∆E) = ρ3 Ẋ 2 τc ∆E0 σ022−∞∼ ρ2 σ0 ∆E0Ẋ .

D ∼(E.23)

Here τc = (σ0 ρẊ)−1 .
From the above arguments, we now know that D ∝ Ẋ 2 for small Ẋ and that D ∝ Ẋfor large Ẋ.

76

References[1] E. Teller, J. Phys. Chem. 41, 109 (1937).
[2] L. D. Landau and E. M. Lifshitz, Quantum Mechanics (Pergamon Press, London,
1958) p. 262.
[3] A. S. Davydov: Quantum Mechanics (Pergamon Press, New York, 1965) p. 175.
[4] J. von Neumann and E. Wigner, Phys. Z. 30, 467. (1929)[5] O. J. Heilmann and E. H. Lieb, Ann. N. Acad. Sci. 172, 583 (1971).
[6] M. V. Berry and M. Wilkinson, Proc. R. Soc. London, Ser. A 392, 15 (1984).
[7] T. Deguchi, K. Fabricius and B. M. McCoy, J. Stat. Phys. 102, 701 (2001).
[8] K. Fabricius and B. M. McCoy, MathPhys Odyssey 2001 ed. M. Kashiwara and T.
Miwa (Birkhäuser, Boston, 2002) p. 119.
[9] T. Deguchi, J. Phys. A 35, 879 (2002).
[10] S. E. Nagler, W. J. L. Buyers, R. L. Armstrong, and B. Briat, 27, 1784 (1983).
[11] H. Kikuchi et al., Physica B 284-288, 1631 (2000).
[12] M. Hagiwara et al., Physica B 294-295, 83 (2001).
[13] T. Deguchi and K. Kudo: Bussei Kenkyu (Kyoto) 80, 195 (2003) [in Japanese].
[14] K. Kudo and T. Deguchi, J. Phys. Soc. Jpn. 72, 1599 (2003).
[15] K. Kudo and T. Deguchi, Phys. Rev. B 68, 052510 (2003).
[16] K. Kudo and T. Deguchi, Phys. Rev. B 69, 132404 (2004).
[17] K. Kudo and T. Deguchi, cond-mat/0409761.
[18] K. Kudo and K. Nakamura, in preparation.

77

[19] A. Klümper, J. R. Reyes Martı́nez, C. Scheeren and M. Shiroishi, J. Stat. Phys. 102,937 (2001).
[20] R. J. Baxter, Exactly Solved Models in Statistical Mechanics (Academic Press, London, 1982).
[21] L. Takhtajan and L. Faddeev, Russ. Math. Survey 34, 11 (1979).
[22] L. A. Takhtajan, Lecture Notes in Physics (Springer-Verlag, Berlin, 1985) Vol. 242,
p. 175.
[23] M. L. Mehta, Random Matrices, 2nd ed. (Academic, San Diego, 1991).
[24] M. V. Berry and M. Tabor, Proc. R. Soc. Lond. A. 356, 375 (1997).
[25] O. Bohigas, M. J. Giannoni, and C. Schmit, Phys. Rev. Lett. 52, 1 (1984).
[26] G. Montambaux, D. Poilblanc, J. Bellissard, and C. Sire, Phys. Rev. Lett. 70, 497(1993).
[27] T. C. Hsu and J. C. Anglès d’Auriac, Phys. Rev. B 47, 14 291 (1993).
[28] D. Poilblanc, T. Ziman, J. Bellissard, F. Mila, and G. Montambaux, Europhys. Lett.
22, 537 (1993).
[29] P. van Ede van der Pals and P. Gaspard, Phys. Rev. E 49, 79 (1994).
[30] J. C. Anglès d’Auriac and J. M. Maillard, Physica A 321, 325 (2003).
[31] D. A. Rabson, B. N. Narozhny, and A. J. Millis, Phys. Rev. B 69, 054403 (2004).
[32] H. Meyer, J. C. Anglès d’Auriac, and H. Bruus, J. Phys. A 29, L483 (1996).
[33] H. Meyer, J. C. Anglès d’Auriac, and J. M. Maillard, Phys. Rev. E 55, 5380 (1997).
[34] H. Meyer and J. C. Anglès d’Auriac, Phys. Rev. E 55, 6608 (1997).
[35] J. C. Anglès d’Auriac and F. Iglói, Phys. Rev. E 58, 241 (1998).
[36] R. Berkovits, and Y. Avishai, J. Phys. C: Condens. Matter 8, 389 (1996).
[37] B. Georgeot and D. L. Shepelyansky, Phys. Rev. Lett. 81, 5129 (1998).
[38] Y. Avishai, J. Richert, and R. Berkovits, Phys. Rev. B 66, 052416 (2002).
[39] L. F. Santos, J. Phys. A 37, 4723 (2004).

78

[40] B. I. Shklovskii, B. Shapiro, B. R. Sears, P. Lambrianides, and H. B. Shore, Phys.
Rev. B 47, 11487 (1993).
[41] N. Nakazono, T. Kato, and K. Nakamura, Prog. Theor. Phys. Supp. 150, 165 (2003);
S. Sawada et al., in preparation.
[42] S. Hirata and K. Nomura, Phys. Rev. B 61, 9453 (2000).
[43] K. Nomura and K. Okamoto, J. Phys. A 27, 5773 (1994).
[44] L. Urba, A. Rosengren, Phys. Rev. B 67, 104406 (2003).
[45] See articles in Quantum Chaos Y2K, Proceedings of Nobel Symposium edited by
K. F. Berggren and S. Åberg (Royal Swedish Academy of Sciences/World Scientiﬁc,
Singapore, 2001).
[46] G. Jona-Lasinio and C. Presilla, Phys. Rev. Lett. 77, 4322 (1996).
[47] T. Prosen, Phys. Rev. Lett. 80, 1808 (1998); T. Prosen, Phys. Rev. E 60, 3949
(1999).
[48] V. V. Flambaum and F. M. Izrailev, Phys. Rev. E 64, 026124 (2001).
[49] M. Wilkinson, J. Phys. A 21, 4021 (1988).
[50] M. Wilkinson and E. J. Austin, Phys. Rev. A 46, 64 (1992).
[51] M. Wilkinson and E. J. Austin, J. Phys. A 28, 2277 (1995).
[52] A. Bulgac, G. D. Dang, and D. Kusnezov, Chaos, Solitons and Fractalts 8, 1149
(1997).
[53] D. Cohen and T. Kottos, Phys. Rev. Lett. 85, 4839 (2000).
[54] D. Cohen, F. M. Izrailev, and T. Kottos, Phys. Rev. Lett. 84, 2052 (2000); D. Cohen,
Lecture Notes in Physics 597 edited by P. Garbaczewski and R. Olkiewicz, p. 317(2002).
[55] M. Machida, K. Saito, and S. Miyashita, J. Phys. Soc. Jpn. 71, 2427 (2002).
[56] Y. Kitaoka et al., J. Phys. Soc. Jpn 67, 3703 (1998); V. P. Plakhty et al., Phys. Rev.
Lett. 85, 3942 (2000); A. Lascialfari et al., Phys. Rev. B 67, 224408 (2003).
[57] K. Nakamura, Y. Nakahara, and A. R. Bishop, Phys. Rev. Lett. 54, 861 (1985);
K. Nakamura, Quantum Chaos (Cambridge University Press, Cambridge, 1993).

79[58] H. Yamasaki, Y. Natsume, A. Terai, and K. Nakamura, J. Phys.: Condens. Matter
16, L395 (2004).
[59] M. Suzuki, Phys. Lett. A 146, 319 (1990).

80