Table Of ContentMatrix product decomposition and classical simulation of quantum dynamics
in the presence of a symmetry
S. Singh,1 H.-Q. Zhou,2 and G. Vidal1
1School of Physical Sciences, the University of Queensland, QLD 4072, Australia
2Department of Physics, Chongqing University, Chongqing 400044, The People’s Republic of China
(Dated: February 4, 2008)
7
0 We propose a refined matrix product state representation for many-body quantum states that
0 are invariant under SU(2) transformations, and indicate how to extend the time-evolving block
2 decimation (TEBD) algorithm in order to simulate time evolution in an SU(2) invariant system.
The resulting algorithm is tested in a critical quantum spin chain and shown to be significantly
n
more efficient than thestandard TEBD.
a
J PACSnumbers:
8
1
Quantummany-bodysystemsaredescribedbyalarge an abelian U(1) symmetry [8]. In addition, it can be
] Hilbert space, one whose dimension grows exponentially cast in the language of spin operators, more familiar to
l
e with the system’s size. This makes the numerical study physicists than group representation theory. As a test,
-
r of generic quantum many-body phenomena computa- we have computed the ground state of the spin-1/2 an-
t
s tionally hard. However, quantum systems are governed tiferromagnetic heisenberg chain, obtaining remarkably
t. by Hamiltonians made of local interactions, that is, by precisetwo-pointcorrelatorsboth forshortandlong dis-
a
highly non-generic operators. As a result, physically rel- tances.
m
evantstates areatypical vectorsin the Hilbert space and In preparation to describe the SU(2) MPS, we start
- may sometimes be described efficiently. Systems in one by introducing a convenient vector basis and discuss a
d
n spatial dimension offer a prominent example. Here the bipartite decomposition of states invariant under SU(2).
o geometry of local interactions induces an anomalously Total spin basis.– Let V be a vector space on which
c small amount of bipartite correlations and an efficient SU(2) acts unitarily by means of transformations ei~v·S~,
[
representation is often possible in terms of a trial wave wherematricesS ,S andS closetheLiealgebrasu(2),
x y z
1 function known as matrix product state (MPS) [1, 2]. namely [S ,S ] = iǫ S , and ~v R3. A total spin
α β αβγ γ
v This, in turn, underlies the success of the density ma- basis (TSB) [V] V satisfies the e∈igenvalue relations
7 trix renormalization group (DMRG) [3], an algorithm to |jtmi∈
2
14 c[4o,m5p,ut6e],girnoculnuddisntgattehs,eatnimdeo-fevsoevlveirnagl rbelcoecnktdeexctiemnasitoionns S~2|[jVtm] i=j(j+1)|[jVtm] i, Sz|[jVtm] i=m|[jVtm] i, (1)
and is associated with the direct sum decomposition of
0 (TEBD) algorithm to simulate time evolution [4].
7 V into irreducible representations (irreps) of SU(2) [9],
Symmetries,offundamentalimportanceinPhysics,re-
0
/ quireaspecialtreatmentinnumericalstudies. Unlessex- V = V˜(j) V(j) . (2)
t plicitlypreservedatthealgorithmiclevel,theyarebound ∼ ⊗
ma to be destroyed by the accumulation of small errors, in Mj (cid:16) (cid:17)
- which case significant features of the system might be Here V˜(j) is a dj-dimensional space that accounts for
d concealed. On the other hand, when properly handled, the degeneracy of the spin-j irrep and has basis [V]
n |jt i ∈
the presence of a symmetry can be exploited to reduce V˜(j), where t = 1, ,d , whereas V(j) is a (2j + 1)-
o j
simulation costs. Whereas the latter has long been re- ···
c dimensionalspace thataccommodatesa spin-j irrepand
v: amliossetdlyinuntheexpcloonrteedxtfoorftDhMe TREGB[3D,7a]l,gtohriethsumbj[e8c].tremains has basis |[jVm]i ∈ V(j), where m is the projection of the
i spin in the z direction, m = j, ,j. Each vector of
X Inthisletterweundertakethestudyofhowtoenhance − ···
the TSB factorizes into degeneracy and irrep parts as
ar theMPSrepresentationandtheTEBDalgorithminsys- [V] = [V] [V] , where Eq. (1) only determines [V] .
tems that are invariant under the action of a Lie group |jtmi |jt i|jmi |jmi
Bipartite decomposition.– A pure state Ψ of a
. We present an explicit theoretical construction of a | i
G bipartite system with vector space A B can always be
refined MPS representation with built-in symmetry, and ⊗
expressed in terms of a TSB for A and a TSB for B as
put forward a significantly faster TEBD algorithm that
both preserves and exploits the symmetry. For simplic- Ψ = Nj1t1m1 [A] [B] . (3)
| i j2t2m2 |j1m1t1i |j2t2m2i
ityandconcreteness,weanalysethesmallestnon-abelian j1Xt1m1j2Xt2m2
case,theSU(2)group,whichisextremelyrelevantinthe
When Ψ is an SU(2) singlet, that is, invariant under
context of isotropic quantum spin systems. The analysis | i
transformations acting simultaneously on A and B, or
of the SU(2) group already contains the major ingredi-
ents of a generic group —in contrast with the case of (S~[A]+S~[B])2 Ψ =0, (S[A]+S[B])Ψ =0, (4)
G | i z z | i
2
thenthesymmetrymaterialisesinconstraintsfortheten- lll ttt ttt
aaa hhh jjj
sor of coefficients N, which splits into degeneracy and aaa aaa jjj ttt jjj
irrep parts according to [10] www jjj
mmm mmm --- mmm
Ψ = Tj [A] [B] ωj [A] [B] , (5) tt tt''
| i Xj Xt1t2 t1t2|jt1i|jt2i! Xm m|jmi|j−mi! aa GG aaaamm'''''' aa '' jj XXjjttjj''tt'' jj''
where ω is completely determined in terms of Clebsch-
Gordan coefficients j j m m j j ;jm [11], namely mm'''' mm ɶɶ mm''
h 1 2 1 2| 1 2 i CCjjmm
mm'''' jj''mm''ssmm''''
(2j+1)−1/2 j =0,1,2,...
ωj , (6)
m ≡ ( 1)m(2j+1)−1/2 j = 1,3,5,...
(cid:26) − 2 2 2 FIG. 1: Diagramatic representation of tensors λ and Γ of an
hj1j2m1m2|j1j2;00i=δj1,j2δm1,−m2ωmj11. (7) MPS and tensors (η,ω) and (X,C˜) of an SU(2) MPS.
Eq. (5)isquitesensible: itsaysthatacoefficientNj1t1m1
j2t2m2 where tensor X relates degeneracy degrees of freedom
in Eq. (3) may be non-zero only if (i) j = j (only the
1 2
and tensor C is given by the Clebsh-Gordan coefficients
product of two spin j irreps can give rise to a spin 0
irrep,that is,the singlet Ψ ) and(ii)m1 = m2,which Cjm = j j m m j j ;jm . (14)
guaranteesthatthez-com|poinentofthespin−vanishes. In j1m1j2m2 h 1 2 1 2| 1 2 i
addition, Eq. (5) embodies the essence of our strategy: MatrixProductdecomposition.–Wenowconsider
toisolatethe degreesoffreedomthatarenotdetermined a chain of n quantum spins with spin s, represented by
by the symmetry – in this case the degeneracy tensor a 1Dlattice where eachsite, labelled by r (r =1,...,n),
Tj . We now considerthe singularvalue decomposition carriesa (2s+1)-dimensionalirrep of SU(2). The coeffi-
t1t2
cients c of a state Ψ of the lattice,
m1m2...mn | i
Tj = (Rj) (ηj) (Sj) (8)
t1t2 t1t t tt2 2s+1 2s+1
Xt |Ψi= ··· cm1m2...mn|[m1]1i|[m2]2i···|[mnn] i, (15)
of tensor Tj for a fixed j, and define mX1=1 mXn=1
t1t2
where [r] is a basis for site r with S[r] [r] = m[r] ,
|Ψ[jAt]i≡ Rt1t|[jAt1]i, |Ψ[jBt]i≡ Stt2|[jBt2]i. (9) can be{c|omdiifi}ed as an MPS [1, 2], z |mi |mi
Xt1 Xt2
c = Γ[1]m1λ[1]Γ[2]m2λ[2] Γ[n]mn. (16)
By combining Eqs. (5), (8) and (9) we arrive to our m1...mn α1 α1 α1α2 α2··· αn−1
canonical symmetric bipartite decomposition (CSBD) α1·X··αn−1
Followingtheconventionsof[4],hereλ[r]aretheSchmidt
α
|Ψi= ηtj|Ψ[jAt]i|Ψ[jBt]i! ωmj |[jAm]i|[jB−]mi!, (10) [cro+effi1ciennts] ooff t|Ψheisapcicnorcdhianign,towthhileebtiepnasrotritΓio[nr]m[1r·e·l·art]es:
Xj Xt Xm ··· αβ
Schmidt vectors for consecutive bipartitions,
which is related to the Schmidt decomposition
2s+1
Ψ = λα Φ[αA] Φ[αB] , (11) |Φ[αr···n]i= Γ[αrβ]mλ[βr] |[mr]i|Φβ[r+1···n]i. (17)
| i | i| i m=1
α X
X
When Ψ is a singlet, that is
by the identifications α (jtm), λ ηjωj and | i
→ α → t m ( S~[r])2 Ψ =0, S[r] Ψ =0, (18)
[A] [A] [A] [B] [B] [B] | i z | i
|Φα i→|Ψjt i|jmi, |Φα i→|Ψjt i|j−mi, (12) Xr Xr
then Eqs. (10) and (13) supersede Eqs. (11) and (17)
where some of the Schmidt coefficients λ are negative.
α
and each tensor λ and Γ in Eq. (16) decomposes into
Moregenerally,astate [CD] ofabipartitesystemC D
|jtm i ⊗ degeneracy and irrep parts, see Fig. (1),
can be expressed in terms of TSBs for C and D as [10]
λ =λ ηj ωj , (19)
α (jtm) → t m
|[jCtmD]i= Xjj1tt1j2t2|j[C1t]1i|[jD2t]2i! Γmαα′′′ =Γ((sjtmm′′))(j′t′m′) → Xjtj′t′ C˜jj′mm′sm′′, (20)
jXj1j2 tXt1t2
where C˜ is related to Clebsch-Gordan coefficients C by
Cjm [C] [D] (13)
mmX1m2 j1m1j2m2|j1m1i|j2m2i! C˜jj′mm′j′′m′′ ≡(−1)2j′(ωmj′′)−1Cjj′mm′j′′m′′. (21)
3
V h[r-1] X[r] h[r] X[r+1] h[r+1] V W
hh [[rr-- 11]] hh [[rr]] hh [[rr++11]] t X1 X2 t' t t'
XX[[rr]] XX[[rr++11]]
j j j j
ww [[rr-- 11]] CC(cid:0)(cid:0)[[rr]] ww [[rr]] CC(cid:0)(cid:0)[[rr++11]] ww [[rr++11]] m VC1 w [r-1]C(cid:1)[r] U[wr,r[+r1]] C(cid:1)[r+1]w [r+1] VC2 - m m w - m
UU[[rr,,rr++11]]
FIG.3: Keystep oftheTEBD algorithm for anSU(2)MPS,
analogous to Figs. (3.i)-(3.iii) in [12] for a regular MPS.—
Once U has been applied on two spins, additional tensors
FIG.2: TheTEBD algorithm isbased onupdatingtheMPS
whenagateU actsontwoneighboringsites. Thisdiagramm VXi,VCi implement a unitary transformation required to re-
absorb these spins into blocks and obtain an updated repre-
generalizesFig. (3.i)in[12]afterthereplacementsλ→(η,ω)
and Γ→(X,C˜) of Eqs. (19)-(20) for an SU(2) MPS. sentation for the bipartition [1···r] : [r+1···n]. Then, for
eachfixedvalueofthej indices(discontinuouslines),theη’s,
X’s and VX’s are multiplied together and the result, with a
weight coming from the product of the ω’s, C˜’s, VC’s and
The SU(2) MPS is defined through Eqs. (19)-(20). In U [that can be pre-computed because none of these tensors
this representation,the constraints imposed by the sym- depend on |Ψi], is added together to give rise to tensor Ω.
metry are used to our advantage. By splitting tensors λ A singular value decomposition of Ωjtjt′ for each value of j
and Γ, we achieve two goals simultaneously. On the one ensues, see Eq. (25), and minor rearrengements finally lead
hand, the resulting MPS is guaranteed, by construction, to updatedtensors X′[r], η′[r] and X′[r+1].
to be invariant under SU(2) transformations. That is,
any algorithm based on this representation will preserve
rithm. We obtain the following comparative costs:
the symmetry exactly and permanently. On the other
hand, all the degrees of freedom of Ψ are concentrated
in smaller tensors η and X (tensors| ωiand C˜ are speci- 3
fiedbythesymmetry),andthustheSU(2)MPSisamore csvd(Θ) (2s+1) [(2j+1)dj] , (24)
∼
economical representation. If denotes the number of j
X
|·|
coefficients of a tensor, then 3
j+2s
csvd(Ω) dj′ . (25)
λ (2j+1)d η d , (22) ∼
| |≡ j → | |≡ j Xj j′≥Xj−2s
j j
X X
Example.– For illustrative purposes, we consider a
|Γ|=(2s+1)|λ||λ′| → |X|≡ djdj′, (23) quantumspinchainwiths=1/2andwithHamiltonian,
(Xj,j′)
H = (S[r]S[r+1]+S[r]S[r+1]+S[r]S[r+1]), (26)
whereλandλ′ arethetensorstotheleftandtotheright x x y y z z
r
ofΓ,andwhere,followingspincompositionrules,thelast X
thatis,thespin-1/2antiferromagneticHeisenbergmodel,
sum is restricted to pairs (j,j′) such that j j′ s.
| − |≤ which is SU(2) invariant and quantum critical at zero
Simulation of time evolution.– Our next step is
temperature. We havecomputedanSU(2)MPSapprox-
to generalize the TEBD algorithm [4] to the simulation
imation to the ground state of H, in the limit n of
of SU(2)-invariant time evolution. This reduces to ex- →∞
aninfinitechain,bysimulatingimaginary-timeevolution
plaining how to update the SU(2) MPS when an SU(2)-
[12] starting from a state made of nearest-neighbor sin-
invariant gate U acts between contiguous sites, see Fig.
glets([r] [r+1] [r] [r+1] )/√2. Withtheconstraint
(2). The update is achieved by following steps analo- |1/2i|−1/2i−|−1/2i|1/2 i
gous to those of the regular TEBD algorithm, see Fig. jdj =600, we have obtained that the following irreps
(3) of [12], involving tensor multiplications and one sin- j, with degeneracies dj, contribute to the odd and even
P
gular value decomposition (SVD), Fig. (3). However, bipartitions [13] of the resulting state,
all these manipulations involve now smaller tensors, and
j 0 1 2 3 4 j 1 3 5 7 9
only tensors X and η of the SU(2) MPS need to be up- 2 2 2 2 2
d 117 247 176 55 5 d 220 242 115 22 1
dated. This results in a substantial reduction of compu- j j
tational space and time, and thus an increase in per-
Eqs. (22)-(25) show substantial computational gains,
formance. For instace, the SVD of Θ in Fig. (3) of
[12], where Θ (2s+1)2 λ2, is now replaced with the Γ 107 csvd(Θ) 9 1010
| |≈ | | | | =50, × =300, (27)
SVD of Ωjtjt′ (see Fig. (3)) for each value of j, where |X| ≈ 2×105 csvd(Ω) ≈ 3×108
|Ωjtjt′| = ( jj+′≥2js−2sdj′)2. The cost csvd(A) of comput- thatis,witharegularMPS,storingthesamestatewould
ing the SVD of a matrix A grows roughly as A3/2 and require about 50 times more computer memory, while
P | |
is the most expensive manipulation of the TEBD algo- performing each SVD would be about 300 times slower.
4
a chain made of large spins, say s = 4, or a spin ladder
e with several legs. We regard the latter as a chain with
n C(r)21100−−87 e13e 2 sSeUv(e2ra)lirsrpeipnssapseirnsEitqe., (w2h)e[r1e0]e.ach site decomposes into
ors i 10−9 ee 4 In addition, the SU(2) MPS is not restricted to the
err 10−10 5e representationofSU(2)singlets. Ontheonehand,itcan
6 be used to represent any SU(2) invariant mixed state ρ
0 1 2 3 4 5 6 7
of the chain, which decomposes as (see Eq. (2))
e c c
10−4 1 1 exact 6c
5
C(r)210−5 e 2 e 3 e 4 e 5 e 6 cc 4 ρ=Mj ρj ⊗I2j+1 (28)
3
10−6 c 2 This is achieved by attaching, to the end of the chain,
an environmentE that duplicates the subspace V of the
500 1,000 2,000 5,000 10,000 20,000
r chain on which ρ is supported, and by considering a sin-
glet purification ΨVE , where ρ = tr ΨVE ΨVE . We
| ρ i E| ρ ih ρ |
firstbuildanSU(2)MPSforthepurificationandthenwe
FIG. 4: (Up) Errors in the two-point correlator C (r) for
2
trace out E. The resulting structure is a matrix product
1 ≤ r ≤ 7, when using an SU(2) MPS of different sizes
χi ∈ {350,700,1110,1450,1800,2200}. Here χ is roughly |λ| representation that retains the advantages of the SU(2)
inEq. (22),thatis,therankofanequivalent(regular) MPS. MPS. In particular, notice that when ρ corresponds to
The lowest line, χ6 =2200, shows theerrors in thedata pre- one single irrep j,
sented in the table. (Down) Numerical results for C (r) for
2
up to r = 20,000 sites, for different sizes χi, together with 1 m
corresponding errors ǫi. ρ= 2j+1 |[jVm]ih[jVm]| (29)
m=−j
X
We have computed the two-point correlators C△(r) then the environment is a site with a spin j, and the
S[0]S[r] and C▽(r) S[1]S[r+1] [14], and the a2verag≡e chain together with the environment is just an extended
h z z i △ 2 ▽≡h z z i spin chain, with the purification being of the form
C (r) (C (r)+C (r))/2 [15]. For small r they read:
2 ≡ 2 2
j
r C2△(r) C2▽(r) C2(r) |Ψρi= √2j1+1 |[jVm]i|[jEm]i. (30)
1 -0.14800224748 -0.14742920605 -0.147715726[7] m=−j
X
2 0.06067976982 0.06067976991 0.060679769[9]
Onthe other hand, the SU(2) MPS canalsobe modified
3 -0.05037860908 -0.05011864581 -0.050248627[4] to represent any pure state [V] of the chain with well
4 0.03465277614 0.03465277645 0.034652776[3] |jmi
defined j and m. To see this, we first consider a mixed
5 -0.0309785296 -0.0308021901 -0.03089036[0] state ρ as in Eq. (29), that is, a symmetrization of [V] ,
6 0.024446726 0.024446726 0.0244467[26] and then a purification Ψ for ρ as in Eq. (29)|,jmfoir
ρ
| i
7 -0.022565932 -0.022430482 -0.0224982[1] whichwecanbuildanSU(2)MPS.Finally,werecallthat
[V] = [E] Ψ , which leads to a simple, SU(2) MPS-
where, for C (r), the square braketsshow the first digits |jmi hjm| ρi
2 [V]
like representation for in terms of the SU(2) MPS
thatdifferfromtheexactsolution[16],fromwhichwere- |jmi
for the purification Ψ . The time-evolution simulation
covere.g. 9significantdigitsforr =1. Anexpressionfor | ρi
techniques described in this paper can be applied to the
thecorrelatorC (2)isalsoknownforlarger[17]. There,
2
above generalized representations.
for r 4,000, 10,000 and 13,000, our results approxi-
≈ Near the completion of this paper, we became aware
mate the asymptotical solution with an error of 1%, 5%
of related results by I. McCulloch derivedindependently
and10%respectively,seeFig. (4). Forcomparison,with
in the context of DMRG [18].
a regular MPS and similar computational resources, we
The authors thank S. Lukyanov for helpful communi-
losethreedigits ofprecisionforr =1,whereasa10%er-
cations. G.V. acknowledgessupport fromthe Australian
rorisalreadyachievedforr 500insteadofr 13,000.
≈ ≈ Research Council through a Federation Fellowship.
Final Remarks.–Theabovetestwithacriticalspin-
1/2 chain unambiguously demonstrates the superiority
of the SU(2) MPS and TEBD with respect to their non-
symmetric versions. Mostpromissingly,these techniques
can now be used to address systems that remain other- [1] M. Fannes, B. Nachtergaele and R. F. Werner, Comm.
wise largely unaccessible to numerical analysis due to a Math. Phys. 144, 3 (1992), pp. 443-490. S. O¨stlund and
largedimensionofthelocalHilbertspace. Theseinclude S. Rommer, Phys. Rev.Lett. 75, 19 (1995), pp.3537.
5
[2] D. Perez-Garcia, F. Verstraete, M.M. Wolf, J.I. Cirac, [10] S. Singh et al, in preparation.
quant-ph/0608197. [11] See,forinstance,http://en.wikipedia.org/wiki/Clebsch-Gordan coefficien
[3] S.R.White,Phys.Rev.Lett.69,2863(1992),Phys.Rev. [12] G. Vidal, cond-mat/0605597.
B48,10345(1993).U.Schollwoeck,Rev.Mod.Phys.77, [13] For semi-integer spins, e.g. s= 1, theCSBD of Eq. (10)
2
259 (2005), cond-mat/0409292. for a partition [1···r]:[r+1···n] of the chain has only
[4] G. Vidal, Phys. Rev. Lett. 91, 147902 (2003), integer (semi-integer) valuesof j when r is even (odd).
quant-ph/0301063.G.Vidal,Phys.Rev.Lett.93,040502 [14] Cw(1) (w =△,▽) are computed by contracting a small
2
(2004), quant-ph/0310089. tensornetworkinvolvingthetensors(η,ω)and(X,C˜)for
[5] A. J. Daley et al, J. Stat. Mech.: Theor. Exp. (2004) twocontiguoussites.Forr>1,Cw(r)iscomputedinthe
2
P04005, cond-mat/0403313. same way, but first we simulate r−1 (SU(2) invariant)
[6] S. R. White and A. E. Feiguin, Phys. Rev. Lett. 93, swap gates that bring thetwo relevant sites together.
076401(2004).F.Verstraete,D.Porras,J.I.CiracPhys. [15] Wefind,numerically,thattheeffectof[13]persistsinthe
Rev.Lett. 93, 227205 (2004); thermodynamic limit (open boundary conditions). As a
[7] S.Rommer,S.Ostlund,Phys.Rev.B55, p.2164 (1997); result, the ground state is invariant under shifts by two
I. P. McCulloch and M. Gul´acsi, Europhys. Lett. 57, (butnotone)latticesites,andC△(r)6=C▽(r)foroddr.
2 2
852 (2002); S. Bergkvist, I. McCulloch, A. Rosengren, [16] See J. Sato, M. Shiroishi and M. Takahashi, Nucl.Phys.
cond-mat/0606265. B 729, 441 (2005), and references therein.
[8] The simple case of an abelian U(1) symmetry was con- [17] Weusec=1.05 inEq.(5.25) ofS.Lukyanov,V.Terras,
sidered by one of us in the context of bridging between Nucl. Phys. B654 (2003) 323-356, hep-th/0206093.
DMRGand TEBD algorithms [5]. [18] I. McCulloch, cond-mat/0701xxx.
[9] N.Landsman, math-ph/9807030.