Improved contraction of finite projected entangled pair states

Basic Information

• Paper ID: 2511.01039
• Title: Improved contraction of finite projected entangled pair states
• Author: Markus Scheb (Ludwig-Maximilians-University Munich)
• Classification: cond-mat.str-el (Strongly Correlated Electron Systems)
• Publication Date: November 4, 2025
• Paper Link: https://arxiv.org/abs/2511.01039

Abstract

This paper presents improved contraction and optimization algorithms for finite projected entangled pair states (fPEPS) and projected entangled pair operators (PEPO). The work comprises two core components: (1) detailed explanation of characteristic contraction patterns appearing in fPEPS calculations and how to slice them to minimize peak memory consumption while ensuring efficient parallel computation; (2) combination of controlled bond expansion (CBE) with randomized singular value decomposition (RSVD), applied throughout the fPEPS algorithm. The paper demonstrates benchmarks on the Hubbard model up to 8×8 system sizes, with SU(2) symmetric bond dimension D=6 (PEPS bonds) and χ=500 (environment bonds).

Research Background and Motivation

Core Problem

Numerical simulation of two-dimensional quantum many-body systems represents a central challenge in condensed matter physics. While one-dimensional systems can be efficiently handled through density matrix renormalization group (DMRG) and matrix product states (MPS), their natural generalization to two dimensions—projected entangled pair states (PEPS)—faces severe computational difficulties.

Problem Significance

1. Fundamental Physics Research: Understanding strongly correlated phenomena such as high-temperature superconductivity and quantum spin liquids requires accurate simulation of two-dimensional quantum systems
2. Methodological Challenges: PEPS suffers from exponentially growing tensor contraction costs due to its loop structure, expectation values cannot be computed exactly, and variational optimization exhibits poor convergence
3. Practical Requirements: Need for unbiased algorithms capable of handling large system sizes, inhomogeneous systems, and open boundary conditions in two-dimensional quantum systems

Limitations of Existing Methods

1. Infinite PEPS (iPEPS): Restricted to small unit cells, cannot handle large inhomogeneous systems
2. Isometric Tensor Networks: While easy to handle, limited to special phases (e.g., string net liquids)
3. Previous fPEPS Implementations: Environment bond dimension χ limited to 250-350, resulting in numerical instability, poor convergence, and large energy errors

Research Motivation

Building upon prior work Ref. 31, the author addresses two main bottlenecks in the fPEPS algorithm:

1. Memory Bottleneck: Intermediate tensor contraction results reach O(χ²(DwD)²), far exceeding other stored tensors
2. Computational Efficiency: The 2-site (2s) algorithm for bond optimization is prohibitively expensive; more efficient 1-site (1s) methods are needed

Core Contributions

1. Optimal Contraction Sequence Design: Detailed exposition of optimal handling of two dominant contraction patterns in fPEPS, using quantum number mapping and slicing techniques to minimize peak memory while ensuring parallel efficiency
2. CBE-RSVD Fusion Framework: First combination of controlled bond expansion (CBE) with randomized singular value decomposition (RSVD), applied to environment approximation and energy minimization, reducing computational cost from 2s to 1s
3. Higher-Fidelity Simulation: Elevation of environment bond dimension from χ=250-350 to χ=500, achieving on the Hubbard model:
   • 4×4 lattice: relative energy error <1% (one supersweep)
   • 6×6 lattice: error reduced from 6.8% to 4.7%
   • 8×8 lattice: error reduced from 16% to 11%, runtime from 21 days to 8 days
4. Practical Validation: Demonstration of stripe structures and antiferromagnetic order in 1/8-filled 8×8 Hubbard model, verifying algorithm's physical reliability

Detailed Methods

Task Definition

Input:

• Two-dimensional lattice Hamiltonian H (represented as PEPO with bond dimension w)
• Initial PEPS wavefunction |ψ⟩ (bond dimension D)
• Environment bond dimension χ

Output:

• Optimized PEPS ground state wavefunction
• Ground state energy E = ⟨ψ|H|ψ⟩

Constraints:

• Open boundary conditions
• Preservation of U(1) or SU(2) symmetry
• Control of computational cost and memory consumption

Core Architecture: fPEPS-PEPO Framework

1. Energy Functional Representation

As shown in Figure 1, the energy functional is represented as a three-layer sandwich structure:

• PEPS Layer: Black bonds connecting adjacent tensors, dimension D
• PEPO Layer: Blue bonds representing the Hamiltonian, dimension w
• Conjugate PEPS Layer: Green bonds connecting physical Hilbert space, dimension d

2. Environment Approximation

Due to exponentially growing exact computation costs from loops, environment approximation is employed: three bond bundles (total dimension DwD) are progressively compressed to cumulative bond dimension χ. In practical simulations, χ ≫ D, w, d.

3. Optimization Scheme

• Local Updates: Scan lattice optimizing individual tensors and adjacent bonds
• Gradient Updates: Simultaneously optimize all PEPS tensors while keeping bond bases fixed
• Supersweeps: Alternating 3 local scans and 100 gradient scans

Technical Innovations

Innovation 1: Optimal Contraction Sequence (Section III)

Problem: Two dominant contraction patterns (Figure 2) directly computed produce enormous intermediate results of size O(χ²(DwD)²)

Solution: Quantum Number Slicing Strategy

1. Preprocessing Phase:
   • Compute and store sandwich tensor C (PEPS-PEPO-PEPS)
   • Scan C, establish quantum number mapping: {(q_t, q_l)}_i → {(q_b, q_r)}_i
   • Build quantum number to dense block mappings for environment tensors T, L, B, R
2. Contraction Execution:
   • Nested loops over (i, q_tr, q_bl)
   • Compute ((T·L)·C)·B (environment approximation) or ((T·L)·C)·(B·R) (effective Hamiltonian action)
   • Contractions for different (i, q_tr, q_bl) are non-overlapping, naturally parallelizable

Advantages:

• Avoids generating enormous intermediate tensors
• Peak memory minimized
• Perfect parallelization without speed loss

Innovation 2: CBE-RSVD Fusion (Section IV)

Background: Bond optimization appears in two stages:

1. Environment approximation (analogous to MPS compression)
2. Energy minimization (analogous to DMRG)

Traditional CBE Problem:

• Requires contraction and decomposition of entire large cluster (Figure 3), cost equivalent to 2s optimization
• "Smart selection" strategy reduces cost but remains inefficient

RSVD Acceleration: For low-rank matrix A of size (χDwD)×(χDwD):

1. Generate (χDwD)×χ̃ Gaussian random matrix Ω (χ̃ ≪ χ)
2. Extract dominant subspace through repeated application of A and A^T to Ω
3. In fPEPS, single AΩ operation suffices

AΩ Operation in Environment Approximation (Figure 5):

• Most expensive contraction cost: O(χ̃χ²(DwD)²)
• Key technique: decompose orthogonal projector into identity and tangential projector
• Avoid explicit computation of tensor leg with dimension χ̄

AΩ Operation in Energy Minimization (Figure 7):

• PEPS tensor processed through weighted trace norm
• χ ≫ D allows direct computation of PEPS orthogonal projector
• Contract Heff|ψ⟩ according to Figure 2(b) structure

Measured Effects:

• 4×4 lattice, SU(2) D=6: reduced from 4 days 150GB (χ=300) to 5 days 26GB (χ=500)
• Memory consumption reduced 83%, simultaneously improving fidelity

Mathematical Details

Completeness Relation (Figure 4): $\mathbb{1} = P_{\text{kept}} + P_{\text{discarded}}$

where χ and χ̄ are dimensions of kept and discarded spaces respectively.

CBE Objective: Under bond dimension χ, maximize overlap of new environment with old environment and PEPS-PEPO-PEPS sandwich, outputting truncated complement space (most important states in discarded space).

Experimental Setup

Test Model: Two-Dimensional Hubbard Model

Hamiltonian: $H = -t\sum_{\langle i,j\rangle,\sigma}(c_{i\sigma}^\dagger c_{j\sigma} + \text{h.c.}) + U\sum_i n_{i\uparrow}n_{i\downarrow}$

Parameter Configuration

• Interaction Strength: U = 8
• Hopping: Nearest-neighbor only
• Boundary Conditions: Open boundaries
• Filling:
  • 4×4 and 6×6: Half-filling (16 and 36 electrons)
  • 8×8: 1/8-filling (56 electrons, inducing stripe structure)

Symmetry Schemes

1. U(1)_spin ⊗ U(1)_charge: D = 4-8
2. SU(2)_spin ⊗ U(1)_charge: D = 4-6

Optimization Strategy

• Supersweeps: 3 local scans + 100 gradient scans
• Environment Bond Dimension: χ = 500 (previous work: 250-400)
• Benchmarks: DMRG calculations (D=4000) energy E_0
  • 4×4: exact ground state energy
  • 6×6 and 8×8: energy upper bounds

Computational Resources

• Runtime: Several days to two weeks
• Memory: 26GB-150GB
• Parallelization: Natural parallelism from quantum number slicing

Experimental Results

Main Results

1. 4×4 Lattice (Figure 8)

Variational Behavior:

• Energy decreases with increasing bond dimension
• Curves flatten as supersweeps increase

Key Improvements (compared to Ref. 31 Figure 38):

• (SU(2), D=4): Now lies between (U(1), D=6) and (U(1), D=7), rather than converging to (U(1), D=5)
• Early Second Supersweep Drop: Indicates that after initial gradient scans, local scan optimization of virtual basis is necessary
• (SU(2), D=6): Breaks 1% error after one supersweep (previous work remained >1% after two supersweeps)

χ Impact:

• χ=500 vs 250-350: Significant convergence improvement
• Low χ not only causes numerical instability but distorts convergence in subtle ways

Performance:

• (SU(2), D=6): 5 days, 26GB (χ=500, 1.5 supersweeps)
• Previous work: 4 days, 150GB (χ=300, 2 supersweeps)

2. 6×6 Lattice (Figure 9)

Variational Behavior: Due to larger system, fewer supersweeps but still showing clear variational properties

Numerical Stability:

• Large bond dimension simulations terminate prematurely (numerical instability)
• Even with computational improvements, rerunning with higher χ remains infeasible
• (U(1), D=7) late-stage temporary increase: temporarily unstable but subsequently recovers

Breakthrough Progress:

• First completion of (U(1), D=8) simulation
• Initially similar to (U(1), D=7), significantly lower in second supersweep

Best Results:

• Previous work: 6.8% (SU(2), D=6)
• This work: 4.7%, higher fidelity (χ=500), runtime halved

3. 8×8 Lattice (Figure 10)

Challenges:

• SU(2) simulation data points sparse (numerical instability or runtime >2 weeks)
• Unlike 4×4 and 6×6, lowest energy comes from U(1) simulation

Best Performance:

• (U(1), D=8): Achieves 11% relative error after one complete supersweep, 8 days
• Previous work: 16% error, half supersweep, 21 days
• Improvement: 31% error reduction, 62% speed improvement

Physical Results Analysis

U(1) Symmetry Local Density (Figure 11)

Stripe Structure (1/8 filling):

• Classical stripe pattern of oscillating charge density
• Incommensurate antiferromagnetic order
• Top/bottom: edge-centered stripes
• Middle: apparently lattice-centered stripes

Limitations: Algorithm far from convergence, cannot determine if this is physical reality or numerical artifact

SU(2) Symmetry Local Density (Figure 12)

• S_z component constructively zero, antiferromagnetic order suppressed
• Shows total spin component
• Charge density exhibits similar stripe structure
• Compared to Ref. 31, distribution more symmetric (qualitative improvement)
• Still requires more supersweeps for reliable physical image

Ablation Studies

While the paper does not explicitly label ablation experiments, different configurations implicitly verify component contributions:

1. χ Value Impact (χ=500 vs 250-350):
   • 4×4: convergence path significantly altered, energy error improved
   • Proves environment bond dimension is critical parameter
2. Symmetry Schemes (U(1) vs SU(2)):
   • Small systems: SU(2) more efficient (same precision, lower D)
   • Large systems: U(1) more stable (8×8)
   • Reveals tradeoff between symmetry and numerical stability
3. Local/Gradient Scan Ratio:
   • Early second supersweep drop suggests 3:100 ratio may be suboptimal
   • Recommends future work revisit this ratio

Related Work

Tensor Network Method Genealogy

One-Dimensional Methods

1. DMRG White 1992-1993: Density matrix renormalization group, gold standard for one-dimensional systems
2. MPS Rommer & Östlund 1997: Tensor network formulation of DMRG, admits canonical form

Two-Dimensional PEPS Variants

1. iPEPS Jordan et al. 2008:
   • Infinite PEPS, successfully applied to various 2D models
   • Limitation: small unit cells
   • Representative work: Corboz's variational optimization, Hubbard model energy extrapolation
2. Isometric Tensor Networks Zaletel & Pollmann 2020:
   • Impose unitarity on virtual bonds
   • Easy to handle but limited to special phases (string net liquids)
3. Gaussian PEPS Emonts & Zohar 2023:
   • Applications to lattice gauge theories, d-wave superconductivity, U(1) Dirac spin liquids
4. PEPS-Monte Carlo Hybrid Vieijra et al. 2021:
   • Direct sampling from PEPS, more efficient contraction but introduces statistical error
   • Liu et al. 2025: exact simulation of Hubbard model
5. Belief Propagation Tindall & Fishman 2023:
   • Attempts to normalize loopy tensor networks
   • Loop series expansion Evenbly et al. 2024: reintroduces neglected contributions

Paper Positioning

• Returns to original fPEPS idea: open boundary conditions without gauge constraints
• Theoretically capable of describing complete physical behavior of large, inhomogeneous 2D quantum systems
• Complements iPEPS (small unit cells) and isometric networks (special gauges)

Technical Inheritance

1. Controlled Bond Expansion (CBE) Gleis et al. 2023: 1s cost bond optimization
2. Randomized SVD Rokhlin et al. 2009; Halko et al. 2011: Fast low-rank matrix decomposition
3. Weighted Trace Norm Evenbly 2018: Construct PEPS orthogonal projector

Conclusions and Discussion

Main Conclusions

1. Technical Breakthroughs:
   • Quantum number slicing strategy achieves memory-optimal, perfectly parallel contraction
   • CBE-RSVD fusion reduces bond optimization cost from 2s to 1s
   • Environment bond dimension elevated from 350 to 500
2. Performance Improvements (compared to Ref. 31):
   • 4×4: energy error <1%, memory reduced 83%
   • 6×6: error reduced from 6.8% to 4.7%, speed doubled
   • 8×8: error reduced from 16% to 11%, speed improved 62%
3. Physical Verification:
   • Successfully reproduced stripe structures in 1/8-filled Hubbard model
   • Observed coexistence of edge-centered and lattice-centered stripes

Limitations

1. Insufficient Competitiveness

Critical Assessment (author's own statement):

"Even with these improvements, energies remain far from DMRG-provided bounds, meaning the fPEPS-PEPO scheme in its current form is still not a competitive tool for computing two-dimensional quantum systems"

Specific Manifestations:

• 6×6 lattice: best error 4.7% (vs DMRG bound)
• 8×8 lattice: 11% error, far from convergence
• Large bond dimension simulations numerically unstable

2. Numerical Stability Issues

• 6×6 and 8×8 large-D simulations frequently terminate prematurely
• SU(2) particularly unstable at 8×8
• Even χ=500 insufficient to guarantee stability

3. Slow Convergence Speed

• 8×8 requires weeks to obtain limited data points
• Local/gradient scan ratio likely suboptimal
• Early second supersweep drop indicates optimization strategy needs adjustment

4. Physical Reliability Questionable

• 8×8 stripe structure may be numerical artifact
• More supersweeps needed for reliable physical conclusions

Future Directions

1. Fusion with Other Methods

Monte Carlo Hybrid Vieijra et al. 2021; Liu et al. 2025:

• More efficient tensor contraction
• Cost: introduces statistical error
• Potential: Liu et al.'s success on Hubbard model

Belief Propagation Tindall & Fishman 2023:

• Attempts to normalize loopy networks
• Does not guarantee convergence to optimal wavefunction
• Loop series expansion Evenbly et al. 2024 can compensate neglected contributions

2. Algorithm Optimization

• Revisit local/gradient scan ratio
• Explore adaptive χ adjustment strategies
• Improve numerical stability (better regularization, preconditioning)

3. Open Questions

Author poses fundamental questions:

"Whether these or other potential new improvements will lead to breakthroughs, and whether fPEPS possesses practical capability to describe arbitrary large, inhomogeneous two-dimensional quantum systems, remain open questions"

In-Depth Evaluation

Strengths

1. Solid Technical Innovation

• Quantum Number Slicing: Elegantly solves memory bottleneck, clear theoretical analysis (Figures 5, 7)
• CBE-RSVD Fusion: Organically combines two mature techniques, significant measured effects
• Complete Implementation Details: Provides sufficient information for others to reproduce

2. Honest Scientific Attitude

• Explicitly states method currently lacks competitiveness
• Objectively analyzes numerical instability
• Maintains caution about physical result reliability
• Rare exemplar of academic integrity

3. Systematic Benchmarking

• Three sizes (4×4, 6×6, 8×8)
• Two symmetries (U(1), SU(2))
• Multiple bond dimensions
• Detailed comparison with prior work, highlighting improvements

4. Clear Visualization

• Numerous tensor network diagrams (Figures 1-7) intuitively display algorithm
• Convergence curves (Figures 8-10) clearly present performance
• Physical results visualized (Figures 11-12)

Shortcomings

1. Fundamental Performance Gap

• Huge gap with DMRG (6×6 error 4.7%, 8×8 error 11%)
• High computational cost (8×8 requires 8 days for only 11% error)
• Fails to demonstrate fPEPS's practical value

2. Numerical Stability Unresolved

• Large-D simulations frequently fail
• χ=500 still insufficient
• Lacks deep analysis of instability origins

3. Experimental Design Limitations

• Lacks Systematic Ablation: Does not separately test slicing strategy and RSVD contributions
• Single Model: Only tests Hubbard model, generalizability unknown
• Insufficient Hyperparameter Exploration: Local/gradient scan ratio not optimized

4. Insufficient Theoretical Analysis

• Why does χ=500 qualitatively improve over 350? Lacks theoretical explanation
• When does numerical instability occur? No predictive indicators
• Is optimal contraction sequence truly optimal? No complexity analysis proof

5. Missing Comparisons with Latest Work

• How does Liu et al. 2025's PEPS-MC perform on Hubbard model?
• Actual effectiveness of belief propagation method?
• Lacks direct comparison with these new methods

Impact Assessment

Contribution to Field

Short-term Impact (Moderate):

• Provides two practical technical improvements for fPEPS community
• Quantum number slicing strategy potentially borrowed by other tensor network methods
• CBE-RSVD fusion framework has general applicability

Long-term Impact (Uncertain):

• Depends on whether fPEPS method ultimately achieves breakthrough
• Author himself skeptical
• May serve more as "how not to do it" lesson

Practical Value

Current (Limited):

• Useful for small-scale research requiring open boundaries, inhomogeneous systems
• 4×4 achieves <1% error, suitable for methodology verification
• 6×6 and above practical utility questionable

Potential (To Be Verified):

• May improve when combined with MC or BP
• As component in more complex hybrid methods

Reproducibility

Strengths:

• Detailed algorithm description, clear diagrams
• Provides specific hyperparameters
• Prior work Ref. 31 provides additional implementation details

Weaknesses:

• Code not open-sourced
• Insufficient symmetry implementation details
• Numerical stability handling not fully explained

Applicable Scenarios

Recommended Scenarios

1. Methodology Research:
   • Test new tensor network techniques
   • Benchmark for more complex methods
   • Small-size (≤4×4) exact verification
2. Specific Physics Problems:
   • Requiring open boundary conditions (iPEPS unsuitable)
   • Inhomogeneous systems (impurities, edge states)
   • Qualitative physics (tolerating few-percent errors)
3. Educational Purposes:
   • Understanding PEPS contraction complexity
   • Learning quantum number symmetry exploitation
   • Tensor network algorithm design case studies

Not Recommended Scenarios

1. High-Precision Requirements: DMRG or QMC more suitable
2. Large-Size Periodic Systems: iPEPS more efficient
3. Production-Level Computing: Current version insufficiently stable and efficient

Methodological Insights

Success Experiences

1. Divide and Conquer: Decompose large problems into manageable subproblems (quantum number slicing)
2. Technical Fusion: Combine existing methods creating synergistic effects (CBE+RSVD)
3. Incremental Improvement: Systematically upgrade upon prior work

Failure Lessons

1. Fundamental Bottlenecks Unbroken: Technical improvements don't address PEPS core difficulties (loops, non-canonicality)
2. Numerical Stability Underestimated: χ's importance only revealed through comparison
3. Performance Evaluation Lagging: Algorithm design insufficiently considered practical application needs

Key References

1. 31 Scheb & Noack, PRB 2023: Direct predecessor, detailed fPEPS-PEPO framework description
2. 32 Gleis et al., PRL 2023: Original controlled bond expansion (CBE) proposal
3. 33-34 Rokhlin et al. 2009; Halko et al. 2011: Theoretical foundations of randomized SVD (RSVD)
4. 36 McCulloch & Osborne 2024: Commentary on CBE, possibly mentions RSVD applications
5. 40 Liu et al., PRL 2025: Latest PEPS-MC hybrid method on Hubbard model
6. 41-42 Tindall & Fishman 2023; Evenbly et al. 2024: Belief propagation and loop series expansion

Summary Commentary

This is a paper of solid technical execution but harsh reality. The author achieves substantial algorithmic engineering improvements (83% memory reduction, 62% speed improvement), but honestly acknowledges these remain insufficient to make fPEPS a practical tool. Such scientific integrity deserves commendation, yet also reveals the profound difficulty of two-dimensional quantum many-body problems.

The paper's value lies more in methodological contributions (quantum number slicing, CBE-RSVD fusion) and sober field assessment rather than solving actual physics problems. It provides useful technical tools and realistic benchmarks for subsequent researchers, but also suggests fPEPS may require deep fusion with other methods (MC, BP) to achieve breakthrough.

For condensed matter theorists, this paper is an excellent tutorial on PEPS method challenges; for algorithm developers, it demonstrates systematic algorithmic improvement; but for physicists needing reliable numerical results, it honestly recommends: currently use DMRG or QMC.