[PPT] - Generalized tensor methods and entanglement measurements for models PowerPoint Presentation

SLIDE 1

Generalized tensor methods and entanglement measurements for models with long-range interactions

¨ Ors Legeza in collaboration with

◮ Jen˝

S´
lyom, Gergely Barcza (Budapest)

◮ Florian Gebhard, Reinhard M. Noack (Marburg) ◮ Valentin Murg, Frank Verstraete (Wien) ◮ Markus Reiher, Konrad H. Marti, Katharina Boguslawski,

Pawel Tecmer, Stefan Knecht (Z¨ urich)

◮ Thorsten R¨

hwedder and Reinhold Schneider (Berlin)

Wigner Research Centre for Physics, Budapest Numerical methods for high-dimensional problems, Paris, 14.04.2014.

SLIDE 2

Brief historical overview

◮ Renormalization Group (Kadanoff transformation) (1966) ◮ Block Renormalization Group (BRG) method ◮ Numerical Renormalization Group (NRG): Wilson(1975) ◮ Density Matrix Renormalization Group (DMRG): White(1992) ◮ Quantum chemistry version of DMRG (QC-DMRG):

White(1999),Mitrushenkov(2001),Chan(2002),Legeza(2002), Reiher(2005),Zgid(2006),Yanai(2008),Xiang(2010), Ma(2012), Wouters(2013)...

◮ Matrix Product State (MPS):

¨ Ostlund, Rommer(1995), Cirac, Verstraete(2004), Oseledets(2009)

SLIDE 3

Some useful reviews, also see additional references therein

◮ U. Schollw¨

ck, Rev. Mod. Phys. 77, 259 (2005).

◮ R. M. Noack and S. R. Manmana, in Diagonalization- and Numerical

Renormalization-Group-Based Methods for Interacting Quantum Systems, (AIP, 2005), vol. 789, pp. 93–163.

◮ F. Verstraete, J.I. Cirac, V. Murg, Adv. Phys. 57 (2), 143 (2008). ◮ ¨

O. Legeza, R. Noack, J. S´
lyom, and L. Tincani, in Computational

Many-Particle Physics, eds. H. Fehske, R. Schneider, and A. Weisse 739, 653–664 (2008).

◮ U. Schollw¨

ck, Ann. Phys. (NY) 326, 96 (2011).

◮ G. K.-L. Chan and S. Sharma, Annu. Rev. Phys. Chem. 62, 465–481

(2011).

◮ K. H. Marti and M. Reiher, Z. Phys. Chem. 224, 583-599 (2010). K. H.

Marti, M. Reiher, Phys. Chem. Chem. Phys. 13 6750-6759 (2011).

◮ G. K.-L. Chan, Wiley Inter-disciplinary Reviews: Computational

Molecular Science 2 (6) (2012) 907–920.

SLIDE 4

Some useful reviews, also see additional references therein

◮ ¨

O. Legeza, T. R¨
hwedder, R. Schneider, Numerical approaches for

high-dimensional PDE’s for quantum chemistry, in Encyclopedia of Applied and Computational Mathematics, Editor-in-chief: Engquist, Bj¨

rn

; Chan, T.; Cook, W.J.; Hairer, E.; Hastad, J.; Iserles, A.; Langtangen, H.P.; Le Bris, C.; Lions, P.L.; Lubich, C.; Majda, A.J.; McLaughlin, J.; Nieminen, R.M.; ODEN, J.; Souganidis, P.; Tveito, A. (Eds.) Springer 2013 ISBN 978-3-540-70530-7

◮ ¨

O. Legeza, T. R¨
hwedder, R. Schneider, Sz. Szalay, arxiv:1310.2736

(2013).

◮ W. Hackbusch, Tensor Spaces and Numerical Tensor Calculus, SSCM

Vol. 42, Springer, 2012.

◮ R. Orus, arXiv:1306.2164. ◮ http://tagung-theoretische-chemie.uni-graz.at/en/past-

workshops/workshop-2014/

SLIDE 5

Topics to be covered

1. Motivation and former approaches:

Model Hamiltonian, i.e., problem to solve Problem in the language of tensor factorization Change of basis and truncation Block Renormalization Group (BRG) Numerical Renormalization Group (NRG)

2. New algorithms in quantum chemistry with polynomial costs:
Density matrix renormalization group (DMRG) White,1992
Matrix Product State (MPS) Ostlund,1995; Verstraete,2004
Tensor Network States (TNS) Marti,2010; Murg,2010; Chan,2013
3. One- and Two-orbital mutual information→Entanglement
Optimizing the algorithms Legeza,2003; Rissler,2006
Efficient construction of active spaces Legeza,2003

Remark: We give a technical introduction to low rank tensor factorization and do not intend to present a detailed review of the field. Only some selected topics will be covered due to time constraint.

SLIDE 6

◮ In ab initio calculations M atoms, Ne electrons,

Coulomb-interaction with = e = m = 1 H =

Ne

i=1

−▽2

i

2 −

Ne,M

i=1,α=1

Zα |ri − rα| + 1 2

Ne

i=j=1

1 |ri − rj| Hamilton-operator, where Zα is the charge of nucleus α, and r denotes the positions of the nuclei and electrons.

◮ The ground state solution is often approximated in a

mean-field manner, where the explicit electron-electron interaction is interchanged by an effective single-particle term. The so called Hartree-Fock (HF) solution of the resulting model often recovers the major part (99 %) of the total energy and provides a good starting point for post-HF computations.

◮ The correlation energy can be defined as the difference of the

exact and the HF energy. The electron correlation, i.e., the beyond mean-field behavior of the electrons, plays an important role in the quantitative description of the chemical properties.

SLIDE 7

Hamiltonian of the interacting electron system

◮ The system is described by the Hamiltonian, for example,

H =

ijσ

Tijσc†

iσcjσ +

ijklσσ′

Vijklσσ′c†

iσc† jσ′ckσ′clσ ◮ Tij denotes the one-electron integral comprising the kinetic

energy of the electrons and the external electric field of the nuclei.

◮ Vijkl stands for the two-electron integrals and contains the e-e

repulsion operator. Vijkl =

d3x1d3x2Φ∗

i (

x1)Φ∗

j (

x2) 1

x1 −

x2 Φk( x2)Φl( x1)

◮ Molecular integrals are calculated via one-electron basis of

atom-centered Gaussians

◮ Major aim: to obtain the desired eigenstates of H.

SLIDE 8

Molecular orbitals, configuration space

Molecular orbitals are obtained, e.g., in a suitable mean-field or MCSCF calculation (for example using the MOLPRO program package). Example.: LiF (6/12)

◮ Two major approximations:

(1) selection of the finite number of basis states (2) restricted number of configurations is taken into account

SLIDE 9

Corrections to the Hartree-Fock state

◮ The full configuration interaction (full-CI) wavefunction can

be expressed in terms of Slater determinants by removing one (S), two (D), three (T) or four (Q) electrons form the HF orb.: ΨFCI = a0ΦSCF +

S

aSΦS +

D

aDΦD +

T

aTΦT + . . . .

◮ Correlation energy: ECorr = EFCI − EHF. ◮ Other expansions are also possible, e.g., Coupled Cluster (see

Schneider’s talk).

SLIDE 10

The problem in the language of tensor factorization

◮ Local tensor space Λ, with dimΛ = q, is denoted by • ◮ Operators for basis |n, σ with n = {0, 1} and σ = {↓, ↑}

c† = 1

, I2 =

1 1

, Ph2 =

1 −1

◮ In a C4 representation • represents a spin-orbital basis with

q = 4. Relevant orbital operators in the |−, | ↓, | ↑, | ↓↑ basis: c†

↓ = c† ⊗I2 =

    1 1     , c†

↑ = Ph2 ⊗c† =

    1 −1     I = I2⊗I2 =     1 1 1 1     , Ph = Ph2⊗Ph2 =     1 −1 −1 1    

SLIDE 11

The problem in the language of tensor factorization

◮ We can put together two C4 tensor spaces, i.e., form a

two-orbital system (••)

◮ Λ(1,2) = Λ(1) ⊗ Λ(2) with dimΛ(1,2) = dimΛ(1)dimΛ(2) = q2 ◮ Basis of the •• system: |φ(1,2) {α1,α2} = |φ(1) α1 ⊗ |φ(2) α2

{α1, α2} = (α1 − 1)q + α2

↓
↑
↓↑

↓ − ↓ ↓ . . . . . . ↓↑ ↓↑ α1, α2 ∈ {1, 2, 3, 4} = {|−, | ↓, | ↑, | ↓↑} {α1α2} ∈ {1, 2, 3, 4 . . . 16}

SLIDE 12

The problem in the language of tensor factorization

◮ Relevant operators for the •• system (dimΛ(1,2) = 4 × 4):

c†

1,↓ = c† ↓ ⊗ I,

c†

2,↓ = Ph ⊗ c† ↓,

I ⇐ I ⊗ I, c†

1,↑ = c† ↑ ⊗ I,

c†

2,↑ = Ph ⊗ c† ↑,

Ph ⇐ Ph ⊗ Ph, H =    h1,1 . . . h1,16 . . . . . . . . . h16,1 . . . h16,16   

◮ Full diagonalization of H gives exact solution (full-CI). ◮ The kth eigenstate of a two-orbital Hamiltonian is

Ψ(1,2)

k

=

α1,α2 U(1,2) α1,α2,k|φ(1) α1 ⊗ |φ(2) α2 ;

α1, α2 = 1 . . . q; k = 1 . . . q2.

SLIDE 13

Change of basis and truncation

◮ We can represent the Hamiltonian in the eigenbasis of the ••

system by a basis state transformation

◮ (Order eigenstates in a descending order, and) form an

perator from the corresponding eigenstates as:

O =       Ψ(1,2)

1

Ψ(1,2)

2.

. . Ψ(1,2)

q2

      , OO† = I

◮ Transform operators to the new basis as H ⇒ OHO†

c†

1,↓ ⇒ Oc† 1,↓O†, c† 2,↓ ⇒ Oc† 2,↓O†, etc ◮ Idea of truncation: select M < q2 states only so OO† = I.

SLIDE 14

Example: two half-spins on ••, |φ↓ ≡ | ↓, |φ↑ ≡ | ↑

Ψm =

α1,α2

Om,2(α1−1)+α2|φα1 ⊗ |φα2, α1, α2 ∈ {1, 2} ≡ {↓, ↑} O ↓ ↓ ↓ ↑ ↑ ↓ ↑ ↑ Sz S Ψ1 1

1

1 Ψ2 1/ √ 2 1/ √ 2 1 Ψ3 1/ √ 2 −1/ √ 2 Ψ4 1 1 1 Take column 1 and 3 for ↓ Take column 2 and 4 for ↑     1 1/ √ 2 −1/ √ 2     = (B(2)[↓])m,α1     1/ √ 2 1/ √ 2 1     = (B(2)[↑])m,α1 In the literature A ≡ BT is used. Ψm =

α ,α

(A(2)[α2])α1,m|φα1 ⊗ |φα2

SLIDE 15

The problem in the language of tensor factorization

◮ For a system with N molecular orbitals: • • • . . . • ◮ Λ(1,2,...,N) = ⊗N i=1Λ(i) with dimΛ(1,2,...,N) = N i dimΛ(i) = qN ◮ Ψ(1,2,...,N) k

=

α1...αN U(1,2,...,N) α1,α2,...αN,k|φ(1) α1 ⊗|φ(2) α2 ⊗. . .⊗|φ(N) αN ; ◮ U(1,2,...,N) α1,α2,...αN,k is a tensor of order N corresponding to the kth

eigenstate of the N-orbital Hamiltonian Example: N = 8

◮ Problem: dimension of U scales exponentially with N

→ We need approximative methods

SLIDE 16

Idea of renormalization: Block Hamiltonians

◮ First we form two-orbital Hamiltonians (blocks) from every

two-orbitals and diagonalize them to obtain their eigenstates: Ψ(1,2)

k

=

α1,α2 U(1,2) α1,α2,k|φ(1) α1 ⊗ |φ(2) α2

Ψ(3,4)

k

=

α3,α4 U(3,4) α3,α4,k|φ(3) α3 ⊗ |φ(4) α4

Ψ(5,6)

k

=

α5,α6 U(5,6) α5,α6,k|φ(5) α5 ⊗ |φ(6) α6

Ψ(7,8)

k

=

α7,α8 U(7,8) α7,α8,k|φ(7) α7 ⊗ |φ(8) α8

α = 1 . . . q and k = 1 . . . q2

SLIDE 17

Idea of renormalization: Truncation

◮ The original idea was that in each RG step we truncate the

Hilbert space of the blocks by keeping only q new states out

f the q2 eigenstates states.

◮ The original form of the Hamiltonian is retained, for lattice

models analytic solution is possible (flow equations, fixed-point, etc (Example: ITF model))

◮ we can keep q < M ≪ qN states but we loose analytic

solution, new operators appear (NRG)

◮ In general we keep dimension of the Hilbert space under

control but we loose information → approximate solution of the problem

◮ Questions:

– how to choose which states to keep at each RG step? – how many states to keep at each RG step? – how accurate will be the truncated solution compared to the full-CI solution?

SLIDE 18

Numerical Renormalization Group (NRG) method

HN−1 HN H2 H H0

0,1

τ

1,2

τ

N−1,N

τ

1

H

2

H

1

Kept states Discarded states

Energy Levels

i i +1 i +2 ...

◮ With exponentially decreasing hopping amplitude, tn = λ−n/2,

where λ is the discretization parameter.

◮ Lowest lying q < M ≪ qN states are kept at each RG

iteration step (works due to separation of energy scales).

◮ Problems for lattice models with λ → 1, blocks are formed

with open boundary condition, segmentation of the total system

SLIDE 19

|m2 =
m1,α2

(A2[α2])α1;m2|α1 ⊗ |α2,

• Use the identity : |m1 =
α1

(A1[α1])1;m1|α1,

• •

|m3 =

m2,α3

(A3[α3])m2;m3;|m2 ⊗ |α3, Transfer tensor: |ml =

ml−1,αl

(Al[αl])ml−1;ml|ml−1⊗|αl, Series of transfer tensors: |ml =

α1,...αl

(A2[α2] . . . Al[αl])α1;ml|α1 . . . αl, For each molecular orbital we can assign a matrix: (Al[αl])ml−1;ml and the wavefunction can be expressed as a product of matrices.

SLIDE 20

Tensor product approximations:

◮ Matrix Product State (MPS) representation:

|Ψ =

q

α1,α2,...,αN

A1

α1A2 α2 · · · AN αN|α1|α2 · · · |αN

Ai[m, m]qi

◮ We can call this a network. Matrix product state (MPS) /

Tensor Train (TT): Ostlund, Rommer (1995), Verstraete, Cirac (2004),

Oseledets (2009), Hackbusch (2009)

SLIDE 21

DMRG provides MPS wavefunction:

Density matrix renormalization group wavefunction: (White, 1992) |ΨTG =

αlαl+1αl+2αr

ψαlαl+1αl+2αr |φ(l)

αl ⊗ |φ(sl) αl+1 ⊗ |φ(sr) αl+2 ⊗ |φ(r) αr

where ψαlαl+1αl+2αr coefficients are determined by an iterative diagonalization of the superblock Hamiltonian. DMRG algorithm provides the optimized set of Ai matrices.

SLIDE 22

SLIDE 23

M

l

Mr q qr

l

B M MR

L

Br

l

B B R

L l

s s r

1. Form and diagonalize the superblock Hamilton operator

|ΨTG =

αlαl+1αl+2αr

ψαlαl+1αl+2αr |φ(l)

αl ⊗|φ(sl) αl+1⊗|φ(sr) αl+2⊗|φ(r) αr

where ψαlαl+1αl+2αr coefficients are determined by an iterative diagonalization of the superblock Hamiltonian.

2. Form a bi-partite representation |ΨTG =

ij ψi,j|φ(L) i

|φ(R)

j

3. Form reduced subsystem density matrix ρ(L)

i,i′ = j ψi,jψ∗ i′,j

4. Diagonalize ρ → ωα eigenvalues, |φ(l)

α eigenstates

5. Form O matrix using M selected |φ(l)

α eigenstates

corresponding to the M largest ωα

6. Renormalize operators: ci ⇒ OciO†

SLIDE 24

Schmidt-decomposition for a bipartite system

◮ For a bipartite system: |ΨT = ij ψij|φL i ⊗ |φR j ◮ Reduced density matrix: ρ(L,R) i,i′

=

j ψijψ∗ i′j ◮ If |Ψ pure state then for |Ψ ∈ Λ = ΛL ⊗ ΛR

|Ψ =

r≤min(ML,MR)

i=1

ωi|ei ⊗ |fi .

◮ |ei, |fi biorthogonal basis, and r is the Schmidt number ◮ If r = 1 ⇒product state, for example, | ↓↑| ↓↑ ◮ If r > 1 ⇒entangled state: non-local property of quantum

mechanics. Example: 1/

√ 2(| ↓| ↑ − | ↑| ↓)

◮ Neumann entropy: s(ργ) = −Tr(ργ ln ργ) , γ ≡ L, R ◮ |ΨT pure state → s(ρL) = s(ρR) ◮ In general, ρL and ρR are in mixed state

SLIDE 25

DMRG wavefunction in MPS form for the l − • • −r superblock Ψ =

{α}
ml
mr

ψmlαl+1αl+2mr ×(Bl[αl] . . . B2[α2])ml;α1 ×(Bl+3[αl+3] . . . BN−1[αN−1])mr;αN × |α1 . . . αN, Connection to the CI-type wavefunction (Marti, Reiher 2011) Ψ =

{α}

C{α}Φ{α}, Therefore, the CI-coefficients in MPS form: C{α} =

Ml

ml

Mr

mr

ψmlαl+1αl+2mr ×(Bl[αl] . . . B2[α2])ml;α1 × (Bl+3[αl+3] . . . BN−1[αN−1])mr;αN

SLIDE 26

Pictorial/diagrammatic description of the one-site DMRG

◮ Component tensors by a dot (or vertex). ◮ Each index or variable by a single line coming out of the vertex ◮ Line connecting two tensors corresponds to an index over

which one has to sum. We call this contraction.

◮ DMRG: on the level of operators; MPS: on the level of states.

5 10 −10 −8 −6 −4 −2 2 4 6 8 10 N or leafs Depth in the tree Tensor−Train

|Ψ〉 〈Ψ|

5 10 −6 −4 −2 2 4 6 N or leafs MPS real site virtual site core tensor 5 10 −6 −4 −2 2 4 6 DMRG Multisite tensorspace

SLIDE 27

Pictorial/diagrammatic description of the one-site DMRG

◮ Component tensors by a dot (or vertex). ◮ Each index or variable by a single line coming out of the vertex ◮ Line connecting two tensors corresponds to an index over

which one has to sum. We call this contraction.

◮ DMRG: on the level of operators; MPS: on the level of states.

5 10 −6 −4 −2 2 4 6 N or leafs Depth in the tree Tensor−Train

|Ψ〉 〈Ψ|

5 10 −6 −4 −2 2 4 6 N or leafs MPS real site virtual site core tensor 5 10 −6 −4 −2 2 4 6 DMRG Multisite tensorspace

SLIDE 28

Pictorial/diagrammatic description of the one-site DMRG

◮ Component tensors by a dot (or vertex). ◮ Each index or variable by a single line coming out of the vertex ◮ Line connecting two tensors corresponds to an index over

which one has to sum. We call this contraction.

◮ DMRG: on the level of operators; MPS: on the level of states.

5 10 −10 −8 −6 −4 −2 2 4 6 8 10 N or leafs Depth in the tree Tensor−Train

|Ψ〉 〈Ψ|

5 10 −6 −4 −2 2 4 6 N or leafs MPS real site virtual site core tensor 5 10 −6 −4 −2 2 4 6 DMRG Multisite tensorspace

SLIDE 29

Error sources and data sparsity

◮ In each DMRG step, the basis states of the system block are

transformed to a new truncated basis set by a unitary transformation based on the preceding SVD

◮ This transformation depends on how accurately the

environment is represented and on the level of truncation

◮ Environmental error, δεsweep, is minimized by a successive

application of the sweepings

◮ Truncation error: δεTR = 1 − M α=1 ωα ◮ For δεsweep → 0, δErel = Const × δεTR (¨ O.L. and G. F´ ath, PRB 1996). ◮ DMRG is a variational method ◮ DMRG is a data-sparse representation of the wavefunction

sparsity ≡ dim(ΛSB)/ dim(ΛFCI)

SLIDE 30

Target state

◮ In the MPS-based approaches, several eigenstates can be

calculated within a single calculation.

◮ Reduced density matrix of the target state, ρ, can be formed

from the reduced density matrices of the lowest n eigenstates as ρ =

γ

pγργ with γ = 1 . . . n,

γ pγ = 1 and Trργ = 1. ◮ Excited states corresponding to the action of given operators

can also be mixed (see Noack’s talk).

SLIDE 31

Example: targeting several states together

5 10 15 20 25 30 35 −106.95 −106.9 −106.85 −106.8 −106.75

E

GS 1XS FCI 5 10 15 20 25 30 35 10

−10

10

−8

10

−6

10

−4

10

∆ Erel

GS 1XS 5 10 15 20 25 30 35 20 40 60

Block states

Ml Mr 5 10 15 20 25 30 35 10

−16

10

−12

10

−8

10

−4

10

δεTr

5 10 15 20 25 30 35 2000 4000 dim(ΛFCI)=134596

dim(ΛSB) Iteration step

SLIDE 32

Dynamic Block State Selection (DBSS) procedure

Optimal truncation scheme: δεTR < ǫ fixed in advance → M is chosen accordingly in every step Example from 2002: DBSS approach applied on F2 (D2h) (18/18) ΛFCI = (2N)!/[(2N − Ne)!Ne!] = 9075135300 ΛDMRG = Ml × 4 × 4 × Mr = 28800000 (Ml = 1200, Mr = 1500) sparsity≃ 315 ΛDMRG = Ml × 4 × 4 × Mr = 7680 (Ml = 120, Mr = 4) sparsity≃ 1181658 ¨

O.L. , Roder, Hess, Phys Rev B (2002)

10 20 30 40 50 60 70 500 1000 1500 2000

Block states F2, L=18, Mmin=256, TREmax=10−8 Ml Mr

10 20 30 40 50 60 70 10

−8

10

−6

10

−4

10

−2

Relative error Ord = [15 10 7 6 16 5 1 3 12 14 13 2 11 4 17 8 18 9]

0.6

SLIDE 33

A priory defined error margin for ground and excited states

One can define the required accuracy prior to the calculations All parameters of the algorithm are adjusted dynamically based on the strength of entanglement encoded in the wavefunction

¨ O.L., S´

lyom, PRB (2003)

SLIDE 34

Part III Entanglement based optimizations

◮ Block entropy →Entanglement

Controlling accuracy
Controlling convergence
Kullback-Leibler relative entropy

◮ One- and Two-orbital mutual information→Entanglement

Optimizing the algorithms
Efficient construction of active spaces
Entanglement and change of basis
Correlation functions and entanglement
Identifying static, dynamic correlations
Description of bond formation and breaking procedures

SLIDE 35

Block entanglement

̺ = |ψ ψ|
̺B = TrA̺

SB = −Tr(̺Bln̺B)

For critical 1 − d systems : sN(l) = c

6 ln 2N π sin πl N

+ g ,

Vidal, Latorre, Rico, Kitaev, PRL (2005), Calabrese, Cardy, JSM (2004), Laflorencie, Sørensen, Chang, Affleck, PRL (2006), ¨ O.L., S´

lyom, Tincani, Noack,

PRL (2007)

SLIDE 36

Mutual information: entanglement correlation

̺ = |ψ ψ|
A subsystem

B subsystem

̺B = TrA̺ SB = −Tr(̺Bln̺B)

p

̺p ⇒ Sp

Sp describes the entanglement of site p with the rest of the system.

p

q

̺p,q ⇒ Sp,q

Sp,q describes the entanglement of site p and q with the rest of the system.

I p,q describes the mutual information between site p and q I p,q = (Sp + Sq − Sp,q)(1 − δp,q)

¨ O.L., S´

lyom, PRB (2003): Quantum Chemistry,

¨ O.L., S´

lyom, PRL (2005): quantum phase transitions (QPT) with q = p + 1.

Rissler, White, Noack, ECP (2005): Quantum chemistry, arbitrary p and q.

SLIDE 37

Block entropy profiles in quantum chemistry (LiF 6/12)

Non-optimized tensor topology Optimized tensor topology

2 4 6 8 10 12 0.1 0.2 Site index Site entropy

01 02 03 04 05

2 4 6 8 10 12 0.1 0.15 0.2 Site index Two−site entropy 2 4 6 8 10 12 0.1 0.2 0.3 Left block length Block entropy 2 3 4 5 6 7 8 9 10 −0.2 −0.1 Left block length Mutual information 2 4 6 8 10 12 0.1 0.2 Site index Site entropy

01 02 03 04 05

2 4 6 8 10 12 0.05 0.1 Site index Two−site entropy 2 4 6 8 10 12 0.05 0.1 0.15 Left block length Block entropy 2 3 4 5 6 7 8 9 10 −0.4 −0.2 Left block length Mutual information

Results as a function of sweepings indicated by different colors
Mutual information: entropy reduction by forming enlarged

blocks:

SLIDE 38

Dynamic Block State Selection (DBSS) procedure

◮ To control the weight of retained information during the RG

procedure: ρ(L) = pkeptρ(L)

kept + (1 − pkept)ρ(L) lost ,

where ρ(L)

kept is formed from the M largest eigenvalues of ρ(L)

and ρ(L)

lost from the remaining eigenvalues with

Trρ(L)

kept = Trρ(L) lost = 1. ◮ The accessible information for such a binary channel would be

less than the Kholevo bound χ ≤ S(ρ) − pkeptS(ρkept) − (1 − pkept)S(ρlost) ,

◮ In DMRG the atypical subspace is neglected → loss of

information.

◮ truncation scheme: χ ≡ S(ρ) − S(ρkept) < ǫ fixed in advance

→ M is chosen accordingly in every step.

¨ O.L., S´

lyom, PRB(2004)

SLIDE 39

◮ Entropy reduction by forming an enlarged system:

sL(l) + sl+1 + IL(l) = sL(l + 1), with IL(l) ≤ 0

◮ Obtained correlation in one RG step:

IL(l) = sL(l + 1) − sL(l) − sl+1 1 ≤ l < N − 1

◮ Entropy sum rule: N−1

l=1

IL(l) = −

N

l=1

sl .

◮ In case of truncation: N−1

l=1

IL(l)+

N

l=1

sl < (N−1)ǫ .

10

−7

10

−6

10

−5

10

−4

10

−3

10

−8

10

−7

10

−6

10

−5

10

−4

10

−3

10

−2

10

−1

χ Relative error CH2 (6/13) H2O (8/24) H2O (10/41) F2 (18/18) SnO (20/59) 10

−7

10

−6

10

−5

10

−4

10

−3

10

−8

10

−7

10

−6

10

−5

10

−4

10

−3

10

−2

10

−1

χ Error in the sum rule

◮ For L = N − 2 → r = ML = MR = 1, correlation functions

can be calculated very quickly

SLIDE 40

Effect of the environment block→efficient warmup

◮ New block states after SVD in each RG step depends on how

accurately the environment is represented.

◮ ρ(N) = ρL ⊗ ρR ◮ For a pure ΨT: sL = sR ◮ Kullback-Leibler relative quantum entropy:

K(ρL||σL) ≡ Tr (ρL ln ρL − ρL ln σL) , where ρL and σL denote reduced density matrices of the left block corresponding to two different right (environment) blocks.

◮ K measures how the system (left) block reduced density

matrix changes as we change the environment.

◮ Optimize environment block based on Kullback-Leibler

entropy.

SLIDE 41

Example: LiF 6/12, orbitals 1,2,3 are the HF orbitals.

Consider three different Mr basis sets:

− − − − − − − − − − − − − − − − − ↓ − − − − − − − − ↑ − − − − − − − − ↓↑ − − − − − − − ↓ − − − − − − − − ↑ − − − − − − − − ↓ ↓ − − − − − − − ↑ ↑ − − − − − − − ↓ ↑ − − − − − − − ↑ ↓ − − − − − − − ↓↑ − − − − − − − ↓ − − − − − − − − ↑ − − − − − − − − ↓ − ↑ − − − − − − ↑ − ↓ − − − − − − ↓ ↓ − − − − − − − ↑ ↑ − − − − − − − − − − − − − − − ↓ − − − − − − − − ↑ − − − − − − − − ↓↑ − − − − − − − − − ↓ − − − − − − − − ↑ − − − − − − − ↓ ↓ − − − − − − − ↑ ↑ − − − − − − − ↑ ↓ − − − − − − − ↓ ↑ − − − − − − − − ↓↑ − − − − − − − − − − ↓ − − − − − − − − ↑ − − − − − ↓ − − ↓ − − − − − ↑ − − ↑ − − − − − ↑ − − ↓ − − − − − ↓ − − ↑ − − − − − − − − − − − − ↓ − − − − − − − − ↑ − − − − − − − − ↓↑ − − − − − − − − − ↓ − − − − − − − − ↑ − − − − − − − ↑ ↓ − − − − − − − ↓ ↑ − − − − − − − − ↓↑ − − − − ↓ − − − − − − − − ↑ − − − − − − − − ↓ − − ↑ − − − − − ↑ − − ↓ − − − − − ↓ − − − ↑ − − − − ↑ − − − ↓ − − − − ↓↑ − − − − − − − −

The reduced density matrix of the L = l + 1 subsystem de- pends on the basis states used for the R = 1+r environment:

2 4 6 8 10 12 14 16 10

−16

10

−14

10

−12

10

−10

10

−8

10

−6

10

−4

10

−2

10

α ωα sL =0.13702 sL =0.21738 sL =4.4515e−05

Case−1 Case−2 Case−3

SLIDE 42

Efficient construction of the active space

The larger the entropy value for a given orbital the larger its contribution to the total correlation energy. Example LiF (6/25):

1 5 10 15 20 25 0.05 0.1 0.15 0.2

Orbital index s(1)

A1 A2 B1 B2

1 5 10 15 20 25 0.2 0.4 0.6 0.8

Orbital index s(1)

A1 A2 B1 B2

r=3.05 r=13.7 CAS-vector ≡ ordering sites with decreasing site entropy values. Include orbitals with largest entropies in the expansion of the active space. CI-based Dynamically Extended Active Space procedure

SLIDE 43

Change of basis and entanglement

H = −t

N

j=1,σ
c+

j,σcj+1,σ + c+ j+1,σcj,σ

+ U

N

j=1

nj,↑nj,↓ H =

kσ

ǫ(k)c†

kσckσ + U

N

k1,k2,k3

c†

k1↑c† k2↓ck3↓ck1+k2−k3↑

ǫ(k) =

r e−ikrt(r), where ki = (2πn)/N, −N/2 < n ≤ N/2.

Tij = −2t cos(ki)δ(i − j) and Vijkl = (U/N)δ(i + j − k − l)

EF k F k F

E

k

2 4 6 8 10 12 14 0.01 0.02 0.03 0.04 0.05

Si U=0.5

2 4 6 8 10 12 14 0.1 0.2

Si U=1

2 4 6 8 10 12 14 0.2 0.4 0.6

Si U=2

2 4 6 8 10 12 14 0.2 0.4 0.6 0.8 1

Si i U=4

real space: U = 0, si = ln 4; U → ∞, si = ln 2. k-space: U = 0, si = 0; U → ∞, si = ln4.

¨ O.L., J. S´

lyom, PRB (2003)

SLIDE 44

Entanglement pattern for a Be ring at r = 2.15 with canonical and localized basis

22 20 24 21 23 19 1 11 6 3 13 8 15 2 16 4 7 12 14 9 17 5 18 10 10 18 5 17 9 14 12 7 4 16 2 15 8 13 3 6 11 1 19 23 21 24 20 22 Orbital index Orbital index 0.005 0.01 0.015 0.02 0.025 0.03 16 11 10 9 17 18 1 2 3 4 5 6 20 21 8 19 13 15 24 23 14 22 12 7 7 12 22 14 23 24 15 13 19 8 21 20 6 5 4 3 2 1 18 17 9 10 11 16 Orbital index Orbital index 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18

22 20 24 21 23 19 1 11 6 3 13 8 15 2 16 4 7 12 14 9 17 5 18 10 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 Orbital index s(1) 16 11 10 9 17 18 1 2 3 4 5 6 20 21 8 19 13 15 24 23 14 22 12 7 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Orbital index s(1)

E. Fertitta, B. Paulus, G. Barcza and ¨

O.L. (2014)

SLIDE 45

Basis states transformation applied to the Hamiltonian

There are two ways to implement the basis transformation: one based on the state and the other based on the Hamiltonian. H =

ij

Tijc†

i cj +

ijkl

Vijklc†

i c† j ckcl ,

The function E(U) can be expressed as E(U) =

ij

˜ T(U)ijc†

i cj +

ijkl

˜ V (U)ijklc†

i c† j ckcl with

˜ T(U) = UTU† ˜ V (U) = (U ⊗ U)V (U ⊗ U)†. The correlation functions c†

i cj and c† i c† j ckcl are calculated with

respect to the original state and are not dependent on the parameters in U. With the function E(U) in this form, its gradient can be calculated explicitly. Both quantities can be evaluated efficiently for different parameter sets U.

Murg, Verstraete, ¨ O.L., Noack (2010)

SLIDE 46

One- (ρi) and two-orbital (ρi,j) reduced density matrix

|ψ =

α1,...,αN

Cα1,...,αN |α1...αN ,

◮ ρi,j is calculated by taking the trace of |ΨΨ| over all local

bases except for αi and αj, the bases of sites i and j, i.e., ρi,j([αi, αj], [α′

i, α′ j]) =

α1,...,✚

✚

αi,...,

✚ ✚

αj,...,αN

Cα1,...,αi,...,αj,...,αN C ∗

α1,...,α′

i,...,α′ j,...,,αN

◮ In the MPS representation, calculation of ρij corresponds to

the contraction of the network except at sites i and j. The computational cost scales as (N − 2)m3q3

◮ This can be decomposed as a sum of projector operators

based on the free variables αi and αj.

◮ ρi and ρi,j can be constructed from operators describing

transitions between single-site basis states.

SLIDE 47

2-site density matrix and generalized correlation functions

Transitions between states of a q-dimensional local Hilbert space: T (m)

i

=

i−1

j=1

I ⊗ T (m) ⊗

L

j=i+1

I. where (T (m))k,l = δ(l+q[k−1]),m for m = 1 . . . q2. Example spin-1/2 boson (qbit): ↓ ↑ ↓ T (1)

i

T (2)

i

↑ T (3)

i

T (4)

i

T (1)

i

−Sz

i + 1 2I

T (2)

i

S−

i

T (3)

i

S+

i

T (4)

i

Sz

i + 1 2I

ρij ↓ ↓ ↓ ↑ ↑ ↓ ↑ ↑ ↓ ↓ T (1)T (1) ↓ ↑ T (1)T (4) T (2)T (3) ↑ ↓ T (3)T (2) T (4)T (1) ↑ ↑ T (4)T (4)

SLIDE 48

Quantum chemistry (some 40 electrons on 40 orbitals)

Example: Task to determine the electronic structure of the binuclear oxo-bridged copper clusters bis(µ-oxo) µ − η2 : η2 peroxo

◮ CASSCF calculations yield a qualitative wrong interpretation

f the energy difference between different isomers.

◮ Too large active space required to get qualitatively correct

picture for standard QC methods.

◮ open d shell problem ◮

K.H. Marti, I.Malkin Ondik, G. Moritz, and M. Reiher, J. Chem. Phys. 128, 014104 (2008).

Y. Kurashige, and T. Yanai, J. Chem. Phys. 130, 234114 (2009).
T. Yanai, Y. Kurashige, E. Neuscamman, and G.K.-L. Chan, J.Chem.Phys. 132, 024105 (2010).
G. Barcza, ¨
O. Legeza, K. H. Marti, and M. Reiher, Phys. Rev. A. 83, 012508 (2011).

SLIDE 49

Site entropy profile→ highly entangled orbitals

bis(µ-oxo) µ − η2 : η2 peroxo

5 10 15 20 25 30 35 40 45 0.2 0.4 0.6 0.8 1 Orbital index s(1) Ag B3u B2u B1g B1u B2g B3g Au 5 10 15 20 25 30 35 40 45 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8

Orbital index s(1)

Ag B3u B2u B1g B1u B2g B3g Au

SLIDE 50

Entanglement picture of the two isomers

bis(µ-oxo) µ − η2 : η2 peroxo

1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 Ag B3u B2u B1g B1u B2g B3g Au 100 10−1 10−2 10−3 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 Ag B3u B2u B1g B1u B2g B3g Au 100 10−1 10−2 10−3

◮ peroxo: orbital pairs 3–14 and 13–35 are highly entangled

→ bonding and anti-bonding orbitals → the O–O bond is intact

◮ bis(µ-oxo): all five orbitals 3, 13, 14, 34, and 35 are entangled

→ four equivalent Cu–O bonds

◮ O–O bond breaking process → transition from the peroxo to the

bisoxo isomer

◮ Mutual information + DMRG → relative energy of the two-isomers

SLIDE 51

Entanglement localization, example for bis(µ-oxo)

energetical ordering

ptimized ordering

1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 Ag B3u B2u B1g B1u B2g B3g Au 100 10−1 10−2 10−3 44 43 42 30 32 11 9 7 21 18 1 5 4 12 16 14 35 36 37 25 24 28 41 40 31 33 29 10 8 20 19 22 2 6 3 17 15 13 34 38 39 26 27 23 Ag B3u B2u B1g B1u B2g B3g Au 100 10−1 10−2 10−3

◮ Reordering orbitals by minimizing the entanglement distance:

ˆ Idist =

i,j Ii,j × |i − j|η , ◮ Apply spectral graph theory: the vector x = (x1, . . . xN) is the

solution that minimizes F(x) = x†Lx =

ij Ii,j(xi − xj)2, with

i xi = 0 and

i x2 i = 1, and the graph Laplacian is

L = D − I with Di,i =

j Ii,j.

The second eigenvector of the Laplacian is the Fiedler vector. Sorting elements of the Fiedler vector → optimal ordering.

SLIDE 52

Entanglement localization

10 20 30 40 0.2 0.4 0.6 0.8 Orbital index Orbital entropy (a) 5 10 15 20 25 30 35 40 0.5 1 Left block length Block entropy (b) 5 10 15 20 25 30 35 40 −1.5 −1 −0.5 Left block length Mutual information (c)

01 02 03 04 05 06

10 20 30 40 0.2 0.4 0.6 0.8 Orbital index Orbital entropy (a) 5 10 15 20 25 30 35 40 0.5 1 Left block length Block entropy (b) 5 10 15 20 25 30 35 40 −1 −0.5 Left block length Mutual information (c)

01 02 03 04 05 06

◮ Height and width of S(l) can be reduced significantly. ◮ Comptational cost is proportional to l S(l).

SLIDE 53

Sorting elements of the Fiedler vector → optimal ordering.

16 11 10 9 17 18 1 2 3 4 5 6 20 21 8 19 13 15 24 23 14 22 12 7 7 12 22 14 23 24 15 13 19 8 21 20 6 5 4 3 2 1 18 17 9 10 11 16 Orbital index Orbital index 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 13 12 1 9 23 14 6 21 24 18 2 19 11 5 16 20 10 4 17 22 7 3 8 15 15 8 3 7 22 17 4 10 20 16 5 11 19 2 18 24 21 6 14 23 9 1 12 13 Orbital index Orbital index 0.02 0.04 0.06 0.08 0.1 0.12

16 11 10 9 17 18 1 2 3 4 5 6 20 21 8 19 13 15 24 23 14 22 12 7 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Orbital index s(1) 13 12 1 9 23 14 6 21 24 18 2 19 11 5 16 20 10 4 17 22 7 3 8 15 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 Orbital index s(1)

Mutual information and site entropy for the ground state of a ring cluster built from six Be atoms at the equilibrium bond length d = 2.15 and at a streched geometry d = 3.30.

SLIDE 54

Relative energy of the isomers: a notebook calculation

50 100 150 200 250 300 −541.54 −541.52 −541.5 −541.48 −541.46 −541.44 −541.42 −541.4 −541.38 −541.36

Energy (a)

Mmin=64, CI=3, opt ORD, opt CAS, Mstart=64 Mmin=64, χ=10−4, CI=3, opt ORD, opt CAS, Mstart=256 Mmin=256, χ=10−5, CI=3, opt ORD, opt CAS, Mstart=512 M=800 in Ref [?]

Fixed number of block states vs DBSS Dynamic Block State Selection sL(l + 1) − sTrunc

L

(l + 1) < χ

◮ DBSS guarantees that the number of block

states are adjusted according to the entanglement between the DMRG blocks and the a priori defined accuracy can be reached.

◮ Dynamically Extended Active Space (CI-DEAS)

method ∆E Reference energies CASSCF(16,14) 1 CASPT2(16,14) 6 bs-B3LYP 221 RASPT2(24,28) 120 Previously published DMRG energies [1],DMRG(26,44)[m=800] 78 [2],DMRG(32,62)[m=2400] 149 [3],DMRG(28,32)[m=2048]-SCF 107 [3],DMRG(28,32)[m=2048]SCF/CT 113 DMRG energies from this work DMRG(26,44)[64/256/10−4] 111 DMRG(26,44)[256/512/10−4] 115 DMRG(26,44)[256/1024/10−4] 113 DMRG(26,44)[256/512/10−5] 113 DMRG(26,44)[256/1024/10−5] 113

M. Reiher et al, J. Chem. Phys. 128, (2008).
T. Yanai et al, J. Chem. Phys. 130, (2009).

G.K.-L. Chan et al, J.Chem.Phys. 132, (2010).

SLIDE 55

Part IV Higher dimensional networks

◮ Tree Tensor Network State (TTNS) algorithm

Multiply connected networks
Structure of the network
Optimization of the network topology
TTNS study of the avoided crossing in LiF

SLIDE 56

Entanglement → Multiply connected networks

LiF at r=3.05 LiF at r=13.7 Rissler, White, Noack, ECP (2006) Murg, Verstraete, Schneider, Nagy, ¨ O.L. (2013) ◮ DMRG → Matrix product states, i.e. optimization along

ne-spatial dimension

◮ Need for an algorithm that reflects the entanglement topology

f the problem → Tensor Network State (TNS) methods

◮ Use tensors Ai[α]m1...mz where z is the coordination number

SLIDE 57

The two-dimensional network

A tree tensor network in which all sites in the tree represent physical orbitals (red lines) and in which entanglement is transferred via the virtual bonds that connect the sites (black lines). Example: each node is represented by a tensor of order zi and the vertical line denotes the physical index αi. The central node is

SLIDE 58

Tree Tensor Network State (TTNS)

|Ψ =

α1,...,αN

Cα1...αN|α1, ..., αN. Cα1...αN describe a tree tensor network, i.e., they emerge from contractions of a set of tensors {A1, . . . , AN}, where Ai [α]m1...mz , is a tensor at each vertex i of the network, with z virtual indices m1 . . . mz of dimension D and one physical index α of dimension q, with z being the coordination number of that site.

A A

(a) (d) (b) (c)

1

m

2

m

3

m

1

m

2

m

3

m

4

m

Vidal, Corboz (2009);

Murg, Verstraete, ¨ O.L., Noack (2010); Nakatani, Chan (2013)

SLIDE 59

Structure of the network

The coefficients Cα1...αN are obtained by contracting the virtual

indices of the tensors.

The structure of the network can be arbitrary.
The coordination number can vary from site to site.
The only condition is that the network is bipartite, i.e., by

cutting one bond, the network separates into two disjoint parts.

For z = 2, the one-dimensional MPS-ansatz used in DMRG is

recovered.

Entanglement is transferred via the virtual bonds that connect

the sites.

For z > 2 the number of virtual bonds required to connect two

arbitrary sites scales logarithmically with the number of sites N, whereas the scaling is linear in N for z = 2.

The maximal distance between two sites, 2∆, scales

logarithmically with N for z > 2.

SLIDE 60

Tree Tensor Network State (TTNS)

site m

r

| h |

1

m,r

h

3 m,r

h

A A A

A
2

m,r

h

(a) (b)

SLIDE 61

Tensor topology optimization:

ij Iij × dη ij

Iij is model dependent:

◮ depends on Tij and Vijkl interaction strengths ◮ depends on the choice of basis ◮ major aim: could we optimize basis on-the-fly (different

approaches are under investigations, Chan, Murg, Verstraete, ¨ O.L., Krumnov, Eiser, Schneider); unsolved problem dij depends on the tensor topology

◮ a possible solution: TTNS with site dependent coordination

number zi.

A A

(a) (d) (b) (c)

1

m

2

m

3

m

1

m

2

m

3

m

4

m

Number of sites in the tree:

N = 1+z

∆

j=1

(z−1)j−1 = z(z − 1)∆ − 2 z − 2 The maximal distance between two orbitals, 2∆, scales loga- rithmically with N for z > 2.

SLIDE 62

Tensor topology optimization:

ij Iij × dη ij (Ex. LiF 6/25)

Energetical ordering (MPS) dij = |i − j| Entanglement localization (MPS) Tree Tensor Network State (TTNS)

SLIDE 63

Optimization of the sweeping

In case of the tree-network, there is more freedom to choose the

ptimal sweeping procedure, i.e., to choose the optimal path

through which the network is traversed. We sweep through the network by going recursively back and forth through each branch. Therefore, according to the labeling of the

rbitals on the lattice shown in the figure one sweep goes through

the orbitals: 1 2 3 4 5 4 6 4 3 7 8 7 3 2 9 10 9 11 9 2 1 12 13 14 13 15 13 12 16 17 16 18 16 12 1 19 20 21 20 22 20 19 23 24 23 25 23 19.

SLIDE 64

TTNS study of the avoided crossing in LiF (6/25)

2 4 6 8 10 12 14 −107.15 −107.1 −107.05 −107 −106.95 −106.9 −106.85 −106.8

Bond length(r) Singlet energies

GS, S=0 1XS, S=0 2XS, S=0 3XS, S=0 2 4 6 8 10 12 14 10

−6

10

−5

10

−4

10

−3

Bond length(r) ∆ E

GS, z=2 GS, z=3 3XS, z=2 3XS, z=3 50 100 150 10

−6

10

−5

10

−4

10

−3

10

−2

Iteration step ∆ EGS

r=3.05 r=5.5 r=11.5 r=13.7 50 100 150 10

−6

10

−5

10

−4

10

−3

10

−2

Iteration step ∆ E3XS

r=3.05 r=5.5 r=11.5 r=13.7 2 4 6 8 10 12 14 1.4 1.6 1.8 2 2.2 2.4 2.6 2.8 3 3.2

Bond length(r) ITot

V. Murg, F. Verstraete, R. Schneider, P. Nagy, ¨
O. L. (2014)