Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
26 commits
Select commit Hold shift + click to select a range
61a88c1
Generalize environment gauge fixing
pbrehmer Jan 21, 2026
0ec910b
Add C4v CTMRG algorithm
pbrehmer Jan 21, 2026
7256f3f
Add C4v CTMRG rrules and defaults
pbrehmer Jan 21, 2026
1a326b9
Fix tensor decompositions
pbrehmer Jan 22, 2026
23c1610
Fix parameter selection
pbrehmer Jan 22, 2026
0a475fe
Make some tests runnable
pbrehmer Jan 22, 2026
094af8c
Fix some more stuff and add fixed_iterscheme C4v test
pbrehmer Jan 22, 2026
a6ce108
Comment out QR things for now
pbrehmer Jan 23, 2026
97760e8
Add gauge fixing and partition function C4v tests
pbrehmer Jan 23, 2026
6effabe
Fix differentiability and add C4v Heisenberg test
pbrehmer Jan 23, 2026
089a005
Automatically symmetrize optimization when alg=:c4v
pbrehmer Jan 23, 2026
f3d89e1
Fix unitcell test imports
pbrehmer Jan 23, 2026
fcdea38
Add eigh_wrapper test and fix IterEigh
pbrehmer Jan 23, 2026
7905570
Fix eigh_wrapper test imports
Yue-Zhengyuan Jan 24, 2026
30f484f
Update src/algorithms/ctmrg/gaugefix.jl
pbrehmer Jan 26, 2026
ee0b5e2
Rename to decomposition algorithm to `decomp_alg`
pbrehmer Jan 26, 2026
f4df969
Fix partition function test import
pbrehmer Jan 26, 2026
0968a4c
Update docstrings to `decomp_alg`
pbrehmer Jan 26, 2026
7651853
Fix flavors test import
pbrehmer Jan 26, 2026
461402e
Add docstrings
pbrehmer Jan 26, 2026
9a7f460
Fix `c4v_projector`
pbrehmer Jan 26, 2026
e2ea8dd
Update docstrings
pbrehmer Jan 28, 2026
3b9d769
Fix `:full` eigh pullback
pbrehmer Jan 28, 2026
2f109f0
Rename `decomp_alg` back to `decomposition_alg`
pbrehmer Jan 30, 2026
e8ea503
Update `initialize_random_c4v_env` and export
pbrehmer Jan 30, 2026
9fb5000
Add singlet corner C4v initialization function
pbrehmer Jan 30, 2026
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion Project.toml
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "PEPSKit"
uuid = "52969e89-939e-4361-9b68-9bc7cde4bdeb"
authors = ["Paul Brehmer", "Lander Burgelman", "Lukas Devos <ldevos98@gmail.com>"]
version = "0.7.0"
authors = ["Paul Brehmer", "Lander Burgelman", "Lukas Devos <ldevos98@gmail.com>"]

[deps]
Accessors = "7d9f7c33-5ae7-4f3b-8dc6-eff91059b697"
Expand Down
11 changes: 11 additions & 0 deletions src/Defaults.jl
Original file line number Diff line number Diff line change
Expand Up @@ -101,9 +101,20 @@ const svd_rrule_alg = :full # ∈ {:full, :gmres, :bicgstab, :arnoldi}
const svd_rrule_broadening = 1.0e-13
const krylovdim_factor = 1.4

# eigh forward & reverse
const eigh_fwd_alg = :qriteration # ∈ {:qriteration, :bisection, :divideandconquer, :multiple, :lanczos, :blocklanczos}
const eigh_rrule_alg = :trunc # ∈ {:trunc, :full}
const eigh_rrule_verbosity = 0

# QR forward & reverse
# const qr_fwd_alg = :something # TODO
# const qr_rrule_alg = :something
# const qr_rrule_verbosity = :something

# Projectors
const projector_alg = :halfinfinite # ∈ {:halfinfinite, :fullinfinite}
const projector_verbosity = 0
const projector_alg_c4v = :c4v_eigh # ∈ {:c4v_eigh, :c4v_qr (TODO)}

# Fixed-point gradient
const gradient_tol = 1.0e-6
Expand Down
9 changes: 7 additions & 2 deletions src/PEPSKit.jl
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,9 @@ include("Defaults.jl") # Include first to allow for docstring interpolation wit

include("utility/util.jl")
include("utility/diffable_threads.jl")
include("utility/eigh.jl")
include("utility/svd.jl")
# include("utility/qr.jl")
include("utility/rotations.jl")
include("utility/hook_pullback.jl")
include("utility/autoopt.jl")
Expand All @@ -42,6 +44,8 @@ include("networks/infinitesquarenetwork.jl")
include("states/infinitepeps.jl")
include("states/infinitepartitionfunction.jl")

include("utility/symmetrization.jl")

include("operators/infinitepepo.jl")
include("operators/transfermatrix.jl")
include("operators/localoperator.jl")
Expand Down Expand Up @@ -71,6 +75,7 @@ include("algorithms/ctmrg/projectors.jl")
include("algorithms/ctmrg/simultaneous.jl")
include("algorithms/ctmrg/sequential.jl")
include("algorithms/ctmrg/gaugefix.jl")
include("algorithms/ctmrg/c4v.jl")

include("algorithms/truncation/truncationschemes.jl")
include("algorithms/truncation/fullenv_truncation.jl")
Expand All @@ -89,8 +94,6 @@ include("algorithms/transfermatrix.jl")
include("algorithms/toolbox.jl")
include("algorithms/correlators.jl")

include("utility/symmetrization.jl")

include("algorithms/optimization/fixed_point_differentiation.jl")
include("algorithms/optimization/peps_optimization.jl")

Expand All @@ -102,6 +105,8 @@ export SVDAdjoint, FullSVDReverseRule, IterSVD
export CTMRGEnv, SequentialCTMRG, SimultaneousCTMRG
export FixedSpaceTruncation, SiteDependentTruncation
export HalfInfiniteProjector, FullInfiniteProjector
export EighAdjoint, IterEigh, C4vCTMRG, C4vEighProjector, C4vQRProjector
export initialize_random_c4v_env, initialize_singlet_c4v_env
export LocalOperator, physicalspace
export product_peps
export reduced_densitymatrix, expectation_value, network_value, cost_function
Expand Down
234 changes: 234 additions & 0 deletions src/algorithms/ctmrg/c4v.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,234 @@
"""
$(TYPEDEF)

CTMRG algorithm assuming a C₄ᵥ-symmetric PEPS, i.e. invariance under 90° spatial rotation and
Hermitian reflection. This requires a single-site unit cell. The projector is obtained from
`eigh` decomposing the Hermitian enlarged corner.

## Fields

$(TYPEDFIELDS)

## Constructors

C4vCTMRG(; kwargs...)

Construct a C₄ᵥ CTMRG algorithm struct based on keyword arguments.
For a full description, see [`leading_boundary`](@ref). The supported keywords are:

* `tol::Real=$(Defaults.ctmrg_tol)`
* `maxiter::Int=$(Defaults.ctmrg_maxiter)`
* `miniter::Int=$(Defaults.ctmrg_miniter)`
* `verbosity::Int=$(Defaults.ctmrg_verbosity)`
* `trunc::Union{TruncationStrategy,NamedTuple}=(; alg::Symbol=:$(Defaults.trunc))`
* `decomposition_alg::Union{<:EighAdjoint,NamedTuple}`
* `projector_alg::Symbol=:$(Defaults.projector_alg_c4v)`
"""
struct C4vCTMRG{P <: ProjectorAlgorithm} <: CTMRGAlgorithm
tol::Float64
maxiter::Int
miniter::Int
verbosity::Int
projector_alg::P
end
function C4vCTMRG(; kwargs...)
return CTMRGAlgorithm(; alg = :c4v, kwargs...)
end
CTMRG_SYMBOLS[:c4v] = C4vCTMRG

"""
$(TYPEDEF)

Projector algorithm implementing the `eigh` decomposition of a Hermitian enlarged corner.

## Fields

$(TYPEDFIELDS)

## Constructors

C4vEighProjector(; kwargs...)

Construct the C₄ᵥ `eigh`-based projector algorithm based on the following keyword arguments:

* `decomposition_alg::Union{<:EighAdjoint,NamedTuple}=EighAdjoint()` : `eigh` algorithm including the reverse rule. See [`EighAdjoint`](@ref).
* `trunc::Union{TruncationStrategy,NamedTuple}=(; alg::Symbol=:$(Defaults.trunc))` : Truncation strategy for the projector computation, which controls the resulting virtual spaces. Here, `alg` can be one of the following:
- `:fixedspace` : Keep virtual spaces fixed during projection
- `:notrunc` : No singular values are truncated and the performed SVDs are exact
- `:truncerror` : Additionally supply error threshold `η`; truncate to the maximal virtual dimension of `η`
- `:truncrank` : Additionally supply truncation dimension `η`; truncate such that the 2-norm of the truncated values is smaller than `η`
- `:truncspace` : Additionally supply truncation space `η`; truncate according to the supplied vector space
- `:trunctol` : Additionally supply singular value cutoff `η`; truncate such that every retained singular value is larger than `η`
* `verbosity::Int=$(Defaults.projector_verbosity)` : Projector output verbosity which can be:
0. Suppress output information
1. Print singular value degeneracy warnings
"""
struct C4vEighProjector{S <: EighAdjoint, T} <: ProjectorAlgorithm
decomposition_alg::S
trunc::T
verbosity::Int
end
function C4vEighProjector(; kwargs...)
return ProjectorAlgorithm(; alg = :c4v_eigh, kwargs...)
end
PROJECTOR_SYMBOLS[:c4v_eigh] = C4vEighProjector

# struct C4vQRProjector{S, T} <: ProjectorAlgorithm
# decomposition_alg::S
# verbosity::Int
# end
# function C4vQRProjector(; kwargs...)
# return ProjectorAlgorithm(; alg = :c4v_qr, kwargs...)
# end
# PROJECTOR_SYMBOLS[:c4v_qr] = C4vQRProjector

#
## C4v-symmetric CTMRG iteration (called through `leading_boundary`)
#

function ctmrg_iteration(
network,
env::CTMRGEnv,
alg::C4vCTMRG,
)
enlarged_corner = c4v_enlarge(network, env, alg.projector_alg)
corner′, projector, info = c4v_projector(enlarged_corner, alg.projector_alg)
edge′ = c4v_renormalize(network, env, projector)
return CTMRGEnv(corner′, edge′), info
end

"""
c4v_enlarge(network, env, ::C4vEighProjector)

Compute the normalized and Hermitian-symmetrized C₄ᵥ enlarged corner.
"""
function c4v_enlarge(network, env, ::C4vEighProjector)
enlarged_corner = TensorMap(EnlargedCorner(network, env, (NORTHWEST, 1, 1)))
return 0.5 * (enlarged_corner + enlarged_corner') / norm(enlarged_corner)
end
# function c4v_enlarge(enlarged_corner, alg::C4vQRProjector)
# # TODO
# end

"""
c4v_projector(enlarged_corner, alg::C4vEighProjector)

Compute the C₄ᵥ projector from `eigh` decomposing the Hermitian `enlarged_corner`.
Also return the normalized eigenvalues as the new corner tensor.
"""
function c4v_projector(enlarged_corner, alg::C4vEighProjector)
trunc = truncation_strategy(alg, enlarged_corner)
D, U, info = eigh_trunc!(enlarged_corner, decomposition_algorithm(alg); trunc)

# Check for degenerate eigenvalues
Zygote.isderiving() && ignore_derivatives() do
if alg.verbosity > 0 && is_degenerate_spectrum(D)
vals = TensorKit.SectorDict(c => diag(b) for (c, b) in blocks(D))
@warn("degenerate eigenvalues detected: ", vals)
end
end

return D / norm(D), U, (; D, U, info...)
end
# function c4v_projector(enlarged_corner, alg::C4vQRProjector)
# # TODO
# end

"""
c4v_renormalize(network, env, projector)

Renormalize the single edge tensor.
"""
function c4v_renormalize(network, env, projector)
new_edge = renormalize_north_edge(env.edges[1], projector, projector', network[1, 1])
return new_edge / norm(new_edge)
end

# TODO: this should eventually be the constructor for a new C4vCTMRGEnv type
function CTMRGEnv(
corner::AbstractTensorMap{T, S, 1, 1}, edge::AbstractTensorMap{T′, S, N, 1}
) where {T, T′, S, N}
corners = fill(corner, 4, 1, 1)
edge_SW = physical_flip(edge)
edges = reshape([edge, edge, edge_SW, edge_SW], (4, 1, 1))
return CTMRGEnv(corners, edges)
end

#
## utility
#

# Adjoint of an edge tensor, but permutes the physical spaces back into the codomain.
# Intuitively, this conjugates a tensor and then reinterprets its 'direction' as an edge tensor.
function _dag(A::AbstractTensorMap{T, S, N, 1}) where {T, S, N}
return permute(A', ((1, (3:(N + 1))...), (2,)))
end
function physical_flip(A::AbstractTensorMap{T, S, N, 1}) where {T, S, N}
return flip(A, 2:N)
end
function project_hermitian(E::AbstractTensorMap{T, S, N, 1}) where {T, S, N}
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can we use a different name for this? I think MatrixAlgebraKit already exports this as well.

E´ = (E + physical_flip(_dag(E))) / 2
return E´
end
function project_hermitian(C::AbstractTensorMap{T, S, 1, 1}) where {T, S}
C´ = (C + C') / 2
return C´
end

# should perform this check at the beginning of `leading_boundary` really...
function check_symmetry(state, ::RotateReflect; atol = 1.0e-10)
@assert length(state) == 1 "check_symmetry only works for single site unit cells"
@assert norm(state[1] - _fit_spaces(rotl90(state[1]), state[1])) /
norm(state[1]) < atol "not rotation invariant"
@assert norm(state[1] - _fit_spaces(herm_depth(state[1]), state[1])) /
norm(state[1]) < atol "not hermitian-reflection invariant"
return nothing
end

#
## environment initialization
#

"""
initialize_random_c4v_env([f=randn, T=scalartype(state)], state, Venv::ElementarySpace)

Initialize a C₄ᵥ-symmetric `CTMRGEnv` on virtual spaces `Venv` with random entries created
by `f` and scalartype `T`.
"""
function initialize_random_c4v_env(state, Venv::ElementarySpace)
return initialize_random_c4v_env(randn, scalartype(state), state, Venv)
end
function initialize_random_c4v_env(f, T, state::InfinitePEPS, Venv::ElementarySpace)
Vpeps = domain(state[1])[1]
return initialize_random_c4v_env(f, T, Vpeps ⊗ Vpeps', Venv)
end
function initialize_random_c4v_env(f, T, state::InfinitePartitionFunction, Venv::ElementarySpace)
Vpf = domain(state[1])[1]
return initialize_random_c4v_env(f, T, Vpf, Venv)
end
function initialize_random_c4v_env(f, T, Vstate::VectorSpace, Venv::ElementarySpace)
corner₀ = DiagonalTensorMap(randn(real(T), Venv ← Venv))
edge₀ = f(T, Venv ⊗ Vstate ← Venv)
edge₀ = project_hermitian(edge₀)
return CTMRGEnv(corner₀, edge₀)
end

"""
initialize_singlet_c4v_env([T=scalartype(state)], state::InfinitePEPS, Venv::ElementarySpace)

Initialize a C₄ᵥ-symmetric `CTMRGEnv` with a singlet corner of dimension `dim(Venv)` and an
identity edge from `id(T, Venv ⊗ Vpeps)`.
"""
function initialize_singlet_c4v_env(state::InfinitePEPS, Venv::ElementarySpace)
return initialize_singlet_c4v_env(scalartype(state), state, Venv)
end
function initialize_singlet_c4v_env(T, state::InfinitePEPS, Venv::ElementarySpace)
Vpeps = domain(state[1])[1]
return initialize_singlet_c4v_env(T, Vpeps, Venv)
end
function initialize_singlet_c4v_env(T, Vpeps::ElementarySpace, Venv::ElementarySpace)
corner₀ = DiagonalTensorMap(zeros(real(T), Venv ← Venv))
corner₀.data[1] = one(real(T))
edge₀ = permute(id(T, Venv ⊗ Vpeps), ((1, 2, 4), (3,)))
return CTMRGEnv(corner₀, edge₀)
end
18 changes: 12 additions & 6 deletions src/algorithms/ctmrg/ctmrg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -19,17 +19,19 @@ function CTMRGAlgorithm(;
maxiter = Defaults.ctmrg_maxiter, miniter = Defaults.ctmrg_miniter,
verbosity = Defaults.ctmrg_verbosity,
trunc = (; alg = Defaults.trunc),
svd_alg = (;),
decomposition_alg = (;),
projector_alg = Defaults.projector_alg, # only allows for Symbol/NamedTuple to expose projector kwargs
)
# replace symbol with projector alg type
haskey(CTMRG_SYMBOLS, alg) || throw(ArgumentError("unknown CTMRG algorithm: $alg"))
alg_type = CTMRG_SYMBOLS[alg]

# parse CTMRG projector algorithm

if alg == :c4v && projector_alg == Defaults.projector_alg
projector_alg = Defaults.projector_alg_c4v
end
projector_algorithm = ProjectorAlgorithm(;
alg = projector_alg, svd_alg, trunc, verbosity
alg = projector_alg, decomposition_alg, trunc, verbosity
)

return alg_type(tol, maxiter, miniter, verbosity, projector_algorithm)
Expand Down Expand Up @@ -64,8 +66,9 @@ supplied via the keyword arguments or directly as an [`CTMRGAlgorithm`](@ref) st
3. Iteration info
4. Debug info
* `alg::Symbol=:$(Defaults.ctmrg_alg)` : Variant of the CTMRG algorithm. See also [`CTMRGAlgorithm`](@ref).
- `:simultaneous`: Simultaneous expansion and renormalization of all sides.
- `:sequential`: Sequential application of left moves and rotations.
- `:simultaneous` : Simultaneous expansion and renormalization of all sides.
- `:sequential` : Sequential application of left moves and rotations.
- `:c4v` : CTMRG assuming C₄ᵥ-symmetric PEPS and environment.

### Projector algorithm

Expand All @@ -76,10 +79,11 @@ supplied via the keyword arguments or directly as an [`CTMRGAlgorithm`](@ref) st
- `:truncrank` : Additionally supply truncation dimension `η`; truncate such that the 2-norm of the truncated values is smaller than `η`
- `:truncspace` : Additionally supply truncation space `η`; truncate according to the supplied vector space
- `:trunctol` : Additionally supply singular value cutoff `η`; truncate such that every retained singular value is larger than `η`
* `svd_alg::Union{<:SVDAdjoint,NamedTuple}` : SVD algorithm for computing projectors. See also [`SVDAdjoint`](@ref). By default, a reverse-rule tolerance of `tol=1e1tol` where the `krylovdim` is adapted to the `env₀` environment dimension.
* `decomposition_alg::Union{<:DecompositionAdjoint,NamedTuple}` : Tensor decomposition algorithm for computing projectors. See e.g. [`SVDAdjoint`](@ref).
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
* `decomposition_alg::Union{<:DecompositionAdjoint,NamedTuple}` : Tensor decomposition algorithm for computing projectors. See e.g. [`SVDAdjoint`](@ref).
* `decomposition_alg` : Tensor decomposition algorithm for computing projectors. See e.g. [`SVDAdjoint`](@ref).

I wonder if it really makes sense to have that type description in the docstring. I don't think most people would know what that type alias means?

* `projector_alg::Symbol=:$(Defaults.projector_alg)` : Variant of the projector algorithm. See also [`ProjectorAlgorithm`](@ref).
- `:halfinfinite` : Projection via SVDs of half-infinite (two enlarged corners) CTMRG environments.
- `:fullinfinite` : Projection via SVDs of full-infinite (all four enlarged corners) CTMRG environments.
- `:c4v_eigh` : Projection via `eigh` of the Hermitian enlarged corner.

## Return values

Expand All @@ -101,6 +105,8 @@ set of vectors and values will be returned as well:
* `U_full` : Last unit cell of all left singular vectors.
* `S_full` : Last unit cell of all singular values.
* `V_full` : Last unit cell of all right singular vectors.

For `C4vCTMRG` instead the last eigendecomposition `U` and `D` (and `U_full`, `D_full`) will be returned.
"""
function leading_boundary(env₀::CTMRGEnv, network::InfiniteSquareNetwork; kwargs...)
alg = select_algorithm(leading_boundary, env₀; kwargs...)
Expand Down
Loading
Loading