diff --git a/src/PEPSKit.jl b/src/PEPSKit.jl index 9bd83c1e7..be6dec4d7 100644 --- a/src/PEPSKit.jl +++ b/src/PEPSKit.jl @@ -114,7 +114,11 @@ include("algorithms/contractions/correlator/pepo_1layer.jl") include("algorithms/ctmrg/sparse_environments.jl") include("algorithms/ctmrg/ctmrg.jl") -include("algorithms/ctmrg/projectors.jl") +include("algorithms/ctmrg/projectors/projectors.jl") +include("algorithms/ctmrg/projectors/halfinfinite.jl") +include("algorithms/ctmrg/projectors/fullinfinite.jl") +include("algorithms/ctmrg/projectors/c4v_eigh.jl") +include("algorithms/ctmrg/projectors/c4v_qr.jl") include("algorithms/ctmrg/simultaneous.jl") include("algorithms/ctmrg/sequential.jl") include("algorithms/ctmrg/gaugefix.jl") diff --git a/src/algorithms/ctmrg/c4v.jl b/src/algorithms/ctmrg/c4v.jl index f8ea37c8b..a4d0da08c 100644 --- a/src/algorithms/ctmrg/c4v.jl +++ b/src/algorithms/ctmrg/c4v.jl @@ -36,75 +36,6 @@ function C4vCTMRG(; kwargs...) end CTMRG_SYMBOLS[:C4vCTMRG] = 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: - - `:FixedSpaceTruncation` : 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 = :C4vEighProjector, kwargs...) -end -PROJECTOR_SYMBOLS[:C4vEighProjector] = C4vEighProjector - -""" -$(TYPEDEF) - -Projector algorithm implementing the `qr` decomposition of a column-enlarged corner. - -## Fields - -$(TYPEDFIELDS) - -## Constructors - - C4vQRProjector(; kwargs...) - -Construct the C₄ᵥ `qr`-based projector algorithm based on the following keyword arguments: - -* `decomposition_alg=QRAdjoint()` : `left_orth` algorithm including the reverse rule. See [`QRAdjoint`](@ref). -""" -struct C4vQRProjector{S} <: ProjectorAlgorithm - # TODO: support all `left_orth` algorithms - decomposition_alg::S -end -function C4vQRProjector(; kwargs...) - return ProjectorAlgorithm(; alg = :C4vQRProjector, kwargs...) -end -PROJECTOR_SYMBOLS[:C4vQRProjector] = C4vQRProjector - -decomposition_algorithm(alg::C4vQRProjector) = alg.decomposition_alg - -# no truncation -_set_truncation(alg::C4vQRProjector, ::TruncationStrategy) = alg -_set_decomposition_truncation(alg::C4vQRProjector, ::TruncationStrategy) = alg - function check_input( ::typeof(leading_boundary), network::InfiniteSquareNetwork, env::CTMRGEnv, alg::C4vCTMRG; atol = 1.0e-10 ) @@ -153,99 +84,22 @@ the corners are identical, and the edges have identical singular values. convergence_tensors(env::CTMRGEnv, ::C4vCTMRG) = (view(env.corners, 1:1, :, :), view(env.edges, 1:1, :, :)) -function ctmrg_iteration(network, env::CTMRGEnv, ::C4vCTMRG{P}) where {P} - throw(ArgumentError("Unknown C4v projector algorithm $P")) -end -function ctmrg_iteration( - network, - env::CTMRGEnv, - alg::C4vCTMRG{<:C4vEighProjector}, - ) - enlarged_corner = c4v_enlarge(network, env, alg.projector_alg) - corner′, projector, info = c4v_projector!(enlarged_corner, alg.projector_alg) +function ctmrg_iteration(network, env::CTMRGEnv, alg::C4vCTMRG) + projector_alg = alg.projector_alg + enlarged_corner = c4v_enlarge(network, env, projector_alg) + projector, info = c4v_projector!(enlarged_corner, projector_alg) edge′ = c4v_renormalize_edge(network, env, projector) - info = (; - contraction_metrics = (; info.truncation_error), - info.D, info.V, - ) - return CTMRGEnv(corner′, edge′), info -end -function ctmrg_iteration( - network, - env::CTMRGEnv, - alg::C4vCTMRG{<:C4vQRProjector}, - ) - enlarged_corner = c4v_enlarge(env, alg.projector_alg) - projector, info = c4v_projector!(enlarged_corner, alg.projector_alg) - edge′ = c4v_renormalize_edge(network, env, projector) - corner′ = c4v_qr_renormalize_corner(edge′, projector, info.R) - info = (; contraction_metrics = (;), info.Q, info.R) + corner′ = c4v_renormalize_corner(edge′, projector, info, projector_alg) 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))) - # TODO: replace by `project_hermitian` - enlarged_corner = 0.5 * (enlarged_corner + enlarged_corner') - return enlarged_corner / norm(enlarged_corner) -end -""" - c4v_enlarge(env, ::C4vQRProjector) - -Compute the normalized column-enlarged northeast corner for C₄ᵥ QR-CTMRG. -""" -function c4v_enlarge(env, ::C4vQRProjector) - return TensorMap(ColumnEnlargedCorner(env, (NORTHWEST, 1, 1))) -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) - alg = _set_decomposition_truncation(alg, truncation_strategy(alg, enlarged_corner)) - eigh_alg = decomposition_algorithm(alg) - - D, V, truncation_error = eigh_trunc!(enlarged_corner, eigh_alg) - - # 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), V, (; D, V, truncation_error) -end -""" - c4v_projector!(enlarged_corner, alg::C4vQRProjector) - -Compute the C₄ᵥ projector by decomposing the column-enlarged corner with `left_orth`. -``` - R--←-- - ↓ - C-←-E-←- = [~Q~] - ↓ | ↓ | -``` -""" -function c4v_projector!(enlarged_corner, alg::C4vQRProjector) - Q, R = left_orth!(enlarged_corner, decomposition_algorithm(alg)) - # TODO: what's a meaningful way to compute a truncation error/condition number in this scheme? - return Q, (; Q, R, truncation_error = zero(scalartype(Q))) +"""Report unsupported projector algorithms when enlarging a C₄ᵥ corner.""" +function c4v_enlarge(network, env, alg::ProjectorAlgorithm) + throw(ArgumentError("Unknown C4v projector algorithm $(typeof(alg))")) end """ - c4v_renormalize_edge(network, env, projector) - -Renormalize the single edge tensor. +Renormalize the single C₄ᵥ edge tensor. ``` |~~~|-←-E-←-|~~~| -←--| P'| | | P |--←- @@ -253,8 +107,8 @@ Renormalize the single edge tensor. | ``` """ -# TODO: possible missing twists for fermions function c4v_renormalize_edge(network, env, projector) + # TODO: possible missing twists for fermions new_edge = renormalize_north_edge(env.edges[1], projector, projector', network[1, 1]) # additional Hermitian projection step for numerical stability new_edge = _project_hermitian(new_edge) @@ -262,25 +116,26 @@ function c4v_renormalize_edge(network, env, projector) end """ - c4v_qr_renormalize_corner(new_edge, projector, R) - -Renormalize the single corner tensor +Renormalize the single C₄ᵥ corner tensor. ``` C-←-E-←-|~~~| | | | P |-←- E---A---|~~~| | | [~P'] - ↓ + ↓ ``` -Using the already calculated QR decomposition + +For `C4vEighProjector`, it is already given by truncated eigenvalues of the enlarged corner. + +For `C4vQRProjector`, we can use the already calculated QR decomposition ``` R--←-- ↓ - C-←-E-←- = [~P~] + C-←-E-←- = [~P~] ↓ | ↓ | ``` -we rewrite the renormalized corner as +to rewrite the renormalized corner as ``` R-←-|~~~| ↓ | P |-←- @@ -290,18 +145,21 @@ we rewrite the renormalized corner as which reuses the renormalized edge `E′` (`new_edge`). (Credit: https://github.com/qiyang-ustc/QRCTM/blob/dd160116c3d7b02076691ceaf0a9833511ae532d/heisenberg.py#L80) """ -# TODO: possible missing twists for fermions -function c4v_qr_renormalize_corner(new_edge::CTMRGEdgeTensor, projector, R) - # contract edge and R +function c4v_renormalize_corner(new_edge::CTMRGEdgeTensor, projector, info, ::C4vEighProjector) + return info.D / norm(info.D) +end +function c4v_renormalize_corner(new_edge::CTMRGEdgeTensor, projector, info, ::C4vQRProjector) + # TODO: possible missing twists for fermions + # contract updated edge and R edge′ = physical_flip(new_edge) - ER = edge′ * twistdual(R, 1) - # contract (edge, R) with projector + ER = edge′ * twistdual(info.R, 1) + # contract (edge′, R) with projector new_corner = contract_edges(ER, projector) new_corner = project_hermitian(new_corner) return new_corner / norm(new_corner) end -# auxilary function: contract two CTMRG edge tensors +"""Contract two CTMRG edge tensors into a corner tensor.""" @generated function contract_edges( EL::CTMRGEdgeTensor{T, S, N}, ER::CTMRGEdgeTensor{T, S, N} ) where {T, S, N} diff --git a/src/algorithms/ctmrg/projectors.jl b/src/algorithms/ctmrg/projectors.jl deleted file mode 100644 index bd81041b1..000000000 --- a/src/algorithms/ctmrg/projectors.jl +++ /dev/null @@ -1,218 +0,0 @@ -""" -$(TYPEDEF) - -Abstract super type for all CTMRG projector algorithms. -""" -abstract type ProjectorAlgorithm end - -const PROJECTOR_SYMBOLS = IdDict{Symbol, Type{<:ProjectorAlgorithm}}() - -""" - ProjectorAlgorithm(; kwargs...) - -Keyword argument parser returning the appropriate `ProjectorAlgorithm` algorithm struct. -""" -function ProjectorAlgorithm(; - alg = Defaults.projector_alg, - decomposition_alg = (;), - trunc = (;), - verbosity = Defaults.projector_verbosity, - ) - # replace symbol with projector alg type - haskey(PROJECTOR_SYMBOLS, alg) || - throw(ArgumentError("unknown projector algorithm: $alg")) - alg_type = PROJECTOR_SYMBOLS[alg] - - # parse SVD forward & rrule algorithm - - decomposition_algorithm = if alg in [:HalfInfiniteProjector, :FullInfiniteProjector] - _alg_or_nt(SVDAdjoint, decomposition_alg) - elseif alg in [:C4vEighProjector] - _alg_or_nt(EighAdjoint, decomposition_alg) - elseif alg in [:C4vQRProjector] - _alg_or_nt(QRAdjoint, decomposition_alg) - end # TODO: how do we solve this in a proper way? - - # qr-ctmrg does not need truncation or degeneracy checks - if alg in [:C4vQRProjector] - return alg_type(decomposition_algorithm) - end - - # parse truncation scheme - truncation_strategy = if trunc isa TruncationStrategy - trunc - elseif trunc isa NamedTuple - _TruncationStrategy(; trunc...) - else - throw(ArgumentError("unknown trunc $trunc")) - end - - return alg_type(decomposition_algorithm, truncation_strategy, verbosity) -end - -""" - decomposition_algorithm(alg::ProjectorAlgorithm) - -Return the tensor decomposition algorithm of the `alg` projector algorithm. -""" -decomposition_algorithm(alg::ProjectorAlgorithm) = alg.decomposition_alg - -function truncation_strategy(alg::ProjectorAlgorithm, edge) - if alg.trunc isa FixedSpaceTruncation - tspace = space(edge, 1) - return isdual(tspace) ? truncspace(flip(tspace)) : truncspace(tspace) - else - return alg.trunc - end -end - -""" - _set_truncation(alg::ProjectorAlgorithm, trunc::TruncationStrategy) - -Update the truncation strategy of a given projector algorithm, keeping all other settings -the same. -""" -function _set_truncation(alg::ProjectorAlgorithm, trunc::TruncationStrategy) - alg = @set alg.trunc = trunc - return alg -end - -""" - _set_decomposition_truncation(alg::ProjectorAlgorithm, trunc::TruncationStrategy) - -Set the truncation strategy of the decomposition algorithm within a given projector -algorithm, keeping all other settings the same. -""" -function _set_decomposition_truncation(alg::ProjectorAlgorithm, trunc::TruncationStrategy) - decomp_alg = decomposition_algorithm(alg) - truncated_fwd_alg = TruncatedAlgorithm(decomp_alg.fwd_alg, trunc) - decomp_alg = @set decomp_alg.fwd_alg = truncated_fwd_alg - alg = @set alg.decomposition_alg = decomp_alg - return alg -end - -""" -$(TYPEDEF) - -Projector algorithm implementing projectors from SVDing the half-infinite CTMRG environment. - -## Fields - -$(TYPEDFIELDS) - -## Constructors - - HalfInfiniteProjector(; kwargs...) - -Construct the half-infinite projector algorithm based on the following keyword arguments: - -* `decomposition_alg::Union{<:SVDAdjoint,NamedTuple}=SVDAdjoint()` : SVD algorithm including the reverse rule. See [`SVDAdjoint`](@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: - - `:FixedSpaceTruncation` : 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 HalfInfiniteProjector{S <: SVDAdjoint, T} <: ProjectorAlgorithm - decomposition_alg::S - trunc::T - verbosity::Int -end -function HalfInfiniteProjector(; kwargs...) - return ProjectorAlgorithm(; alg = :HalfInfiniteProjector, kwargs...) -end - -PROJECTOR_SYMBOLS[:HalfInfiniteProjector] = HalfInfiniteProjector - -""" -$(TYPEDEF) - -Projector algorithm implementing projectors from SVDing the full 4x4 CTMRG environment. - -## Fields - -$(TYPEDFIELDS) - -## Constructors - - FullInfiniteProjector(; kwargs...) - -Construct the full-infinite projector algorithm based on the following keyword arguments: - -* `decomposition_alg::Union{<:SVDAdjoint,NamedTuple}=SVDAdjoint()` : SVD algorithm including the reverse rule. See [`SVDAdjoint`](@ref). -* `trunc::Union{TruncationStrategy,NamedTuple}=(; alg::Symbol=:$(Defaults.trunc))` : Truncation scheme for the projector computation, which controls the resulting virtual spaces. Here, `alg` can be one of the following: - - `:FixedSpaceTruncation` : 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 FullInfiniteProjector{S <: SVDAdjoint, T} <: ProjectorAlgorithm - decomposition_alg::S - trunc::T - verbosity::Int -end -function FullInfiniteProjector(; kwargs...) - return ProjectorAlgorithm(; alg = :FullInfiniteProjector, kwargs...) -end - -PROJECTOR_SYMBOLS[:FullInfiniteProjector] = FullInfiniteProjector - -""" - compute_projector(enlarged_corners, alg::ProjectorAlgorithm) - -Determine left and right projectors at the bond given determined by the enlarged corners -using the specified `alg`. -""" -function compute_projector(enlarged_corners, alg::HalfInfiniteProjector) - # SVD half-infinite environment - halfinf = half_infinite_environment(enlarged_corners...) - svd_alg = decomposition_algorithm(alg) - U, S, V, truncation_error = svd_trunc!(halfinf / norm(halfinf), svd_alg) - - # get some decomposition info - truncation_error = truncation_error / norm(S) # normalize truncation error - - # Check for degenerate singular values - Zygote.isderiving() && ignore_derivatives() do - if alg.verbosity > 0 && is_degenerate_spectrum(S) - svals = TensorKit.SectorDict(c => diag(b) for (c, b) in blocks(S)) - @warn("degenerate singular values detected: ", svals) - end - end - - P_left, P_right = contract_projectors(U, S, V, enlarged_corners...) - return (P_left, P_right), (; U, S, V, truncation_error) -end -function compute_projector(enlarged_corners, alg::FullInfiniteProjector) - halfinf_left = half_infinite_environment(enlarged_corners[1], enlarged_corners[2]) - halfinf_right = half_infinite_environment(enlarged_corners[3], enlarged_corners[4]) - - # SVD full-infinite environment - fullinf = full_infinite_environment(halfinf_left, halfinf_right) - svd_alg = decomposition_algorithm(alg) - U, S, V, truncation_error = svd_trunc!(fullinf / norm(fullinf), svd_alg) - - # get some decomposition info - truncation_error = truncation_error / norm(S) # normalize truncation error - - # Check for degenerate singular values - Zygote.isderiving() && ignore_derivatives() do - if alg.verbosity > 0 && is_degenerate_spectrum(S) - svals = TensorKit.SectorDict(c => diag(b) for (c, b) in blocks(S)) - @warn("degenerate singular values detected: ", svals) - end - end - - P_left, P_right = contract_projectors(U, S, V, halfinf_left, halfinf_right) - return (P_left, P_right), (; U, S, V, truncation_error) -end diff --git a/src/algorithms/ctmrg/projectors/c4v_eigh.jl b/src/algorithms/ctmrg/projectors/c4v_eigh.jl new file mode 100644 index 000000000..9e4903900 --- /dev/null +++ b/src/algorithms/ctmrg/projectors/c4v_eigh.jl @@ -0,0 +1,72 @@ +""" +$(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: + - `:FixedSpaceTruncation` : 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 = :C4vEighProjector, kwargs...) +end +PROJECTOR_SYMBOLS[:C4vEighProjector] = C4vEighProjector + +""" +Compute the normalized and Hermitian-symmetrized C₄ᵥ enlarged corner. +``` + C-←-E-←- + | | + E---A--- + | | +``` +""" +function c4v_enlarge(network, env, ::C4vEighProjector) + enlarged_corner = TensorMap(EnlargedCorner(network, env, (NORTHWEST, 1, 1))) + enlarged_corner = project_hermitian(enlarged_corner) + return enlarged_corner / norm(enlarged_corner) +end + +""" +Compute the C₄ᵥ projector from `eigh` decomposing the Hermitian `enlarged_corner`. +Return the projector and decomposition diagnostics used to renormalize the corner. +""" +function c4v_projector!(enlarged_corner, alg::C4vEighProjector) + alg = _set_decomposition_truncation(alg, truncation_strategy(alg, enlarged_corner)) + eigh_alg = decomposition_algorithm(alg) + + D, V, truncation_error = eigh_trunc!(enlarged_corner, eigh_alg) + + # 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 V, (; contraction_metrics = (; truncation_error), D, V) +end diff --git a/src/algorithms/ctmrg/projectors/c4v_qr.jl b/src/algorithms/ctmrg/projectors/c4v_qr.jl new file mode 100644 index 000000000..51dfbd4f8 --- /dev/null +++ b/src/algorithms/ctmrg/projectors/c4v_qr.jl @@ -0,0 +1,57 @@ +""" +$(TYPEDEF) + +Projector algorithm implementing the `qr` decomposition of a column-enlarged corner. + +## Fields + +$(TYPEDFIELDS) + +## Constructors + + C4vQRProjector(; kwargs...) + +Construct the C₄ᵥ `qr`-based projector algorithm based on the following keyword arguments: + +* `decomposition_alg=QRAdjoint()` : `left_orth` algorithm including the reverse rule. See [`QRAdjoint`](@ref). +""" +struct C4vQRProjector{S} <: ProjectorAlgorithm + # TODO: support all `left_orth` algorithms + decomposition_alg::S +end +function C4vQRProjector(; kwargs...) + return ProjectorAlgorithm(; alg = :C4vQRProjector, kwargs...) +end +PROJECTOR_SYMBOLS[:C4vQRProjector] = C4vQRProjector + +decomposition_algorithm(alg::C4vQRProjector) = alg.decomposition_alg + +# no truncation +_set_truncation(alg::C4vQRProjector, ::TruncationStrategy) = alg +_set_decomposition_truncation(alg::C4vQRProjector, ::TruncationStrategy) = alg + +""" +Compute the column-enlarged northwest corner for C₄ᵥ QR-CTMRG. +``` + C-←-E-←- + ↓ | +``` +""" +function c4v_enlarge(network, env, ::C4vQRProjector) + return TensorMap(ColumnEnlargedCorner(env, (NORTHWEST, 1, 1))) +end + +""" +Compute the C₄ᵥ projector by decomposing the column-enlarged corner with `left_orth`. +``` + R--←-- + ↓ + C-←-E-←- = [~Q~] + ↓ | ↓ | +``` +""" +function c4v_projector!(enlarged_corner, alg::C4vQRProjector) + Q, R = left_orth!(enlarged_corner, decomposition_algorithm(alg)) + # TODO: what's a meaningful way to compute a truncation error/condition number in this scheme? + return Q, (; contraction_metrics = (;), Q, R) +end diff --git a/src/algorithms/ctmrg/projectors/fullinfinite.jl b/src/algorithms/ctmrg/projectors/fullinfinite.jl new file mode 100644 index 000000000..61e7dc3cb --- /dev/null +++ b/src/algorithms/ctmrg/projectors/fullinfinite.jl @@ -0,0 +1,61 @@ +""" +$(TYPEDEF) + +Projector algorithm implementing projectors from SVDing the full 4x4 CTMRG environment. + +## Fields + +$(TYPEDFIELDS) + +## Constructors + + FullInfiniteProjector(; kwargs...) + +Construct the full-infinite projector algorithm based on the following keyword arguments: + +* `decomposition_alg::Union{<:SVDAdjoint,NamedTuple}=SVDAdjoint()` : SVD algorithm including the reverse rule. See [`SVDAdjoint`](@ref). +* `trunc::Union{TruncationStrategy,NamedTuple}=(; alg::Symbol=:$(Defaults.trunc))` : Truncation scheme for the projector computation, which controls the resulting virtual spaces. Here, `alg` can be one of the following: + - `:FixedSpaceTruncation` : 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 FullInfiniteProjector{S <: SVDAdjoint, T} <: ProjectorAlgorithm + decomposition_alg::S + trunc::T + verbosity::Int +end +function FullInfiniteProjector(; kwargs...) + return ProjectorAlgorithm(; alg = :FullInfiniteProjector, kwargs...) +end + +PROJECTOR_SYMBOLS[:FullInfiniteProjector] = FullInfiniteProjector + +function compute_projector(enlarged_corners, alg::FullInfiniteProjector) + halfinf_left = half_infinite_environment(enlarged_corners[1], enlarged_corners[2]) + halfinf_right = half_infinite_environment(enlarged_corners[3], enlarged_corners[4]) + + # SVD full-infinite environment + fullinf = full_infinite_environment(halfinf_left, halfinf_right) + svd_alg = decomposition_algorithm(alg) + U, S, V, truncation_error = svd_trunc!(fullinf / norm(fullinf), svd_alg) + + # get some decomposition info + truncation_error = truncation_error / norm(S) # normalize truncation error + + # Check for degenerate singular values + Zygote.isderiving() && ignore_derivatives() do + if alg.verbosity > 0 && is_degenerate_spectrum(S) + svals = TensorKit.SectorDict(c => diag(b) for (c, b) in blocks(S)) + @warn("degenerate singular values detected: ", svals) + end + end + + P_left, P_right = contract_projectors(U, S, V, halfinf_left, halfinf_right) + return (P_left, P_right), (; U, S, V, truncation_error) +end diff --git a/src/algorithms/ctmrg/projectors/halfinfinite.jl b/src/algorithms/ctmrg/projectors/halfinfinite.jl new file mode 100644 index 000000000..9b79c8260 --- /dev/null +++ b/src/algorithms/ctmrg/projectors/halfinfinite.jl @@ -0,0 +1,64 @@ +""" +$(TYPEDEF) + +Projector algorithm implementing projectors from SVDing the half-infinite CTMRG environment. + +## Fields + +$(TYPEDFIELDS) + +## Constructors + + HalfInfiniteProjector(; kwargs...) + +Construct the half-infinite projector algorithm based on the following keyword arguments: + +* `decomposition_alg::Union{<:SVDAdjoint,NamedTuple}=SVDAdjoint()` : SVD algorithm including the reverse rule. See [`SVDAdjoint`](@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: + - `:FixedSpaceTruncation` : 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 HalfInfiniteProjector{S <: SVDAdjoint, T} <: ProjectorAlgorithm + decomposition_alg::S + trunc::T + verbosity::Int +end +function HalfInfiniteProjector(; kwargs...) + return ProjectorAlgorithm(; alg = :HalfInfiniteProjector, kwargs...) +end + +PROJECTOR_SYMBOLS[:HalfInfiniteProjector] = HalfInfiniteProjector + +""" + compute_projector(enlarged_corners, alg::ProjectorAlgorithm) + +Determine left and right projectors at the bond given determined by the enlarged corners +using the specified `alg`. +""" +function compute_projector(enlarged_corners, alg::HalfInfiniteProjector) + # SVD half-infinite environment + halfinf = half_infinite_environment(enlarged_corners...) + svd_alg = decomposition_algorithm(alg) + U, S, V, truncation_error = svd_trunc!(halfinf / norm(halfinf), svd_alg) + + # get some decomposition info + truncation_error = truncation_error / norm(S) # normalize truncation error + + # Check for degenerate singular values + Zygote.isderiving() && ignore_derivatives() do + if alg.verbosity > 0 && is_degenerate_spectrum(S) + svals = TensorKit.SectorDict(c => diag(b) for (c, b) in blocks(S)) + @warn("degenerate singular values detected: ", svals) + end + end + + P_left, P_right = contract_projectors(U, S, V, enlarged_corners...) + return (P_left, P_right), (; U, S, V, truncation_error) +end diff --git a/src/algorithms/ctmrg/projectors/projectors.jl b/src/algorithms/ctmrg/projectors/projectors.jl new file mode 100644 index 000000000..6815669f8 --- /dev/null +++ b/src/algorithms/ctmrg/projectors/projectors.jl @@ -0,0 +1,92 @@ +""" +$(TYPEDEF) + +Abstract super type for all CTMRG projector algorithms. +""" +abstract type ProjectorAlgorithm end + +const PROJECTOR_SYMBOLS = IdDict{Symbol, Type{<:ProjectorAlgorithm}}() + +""" + ProjectorAlgorithm(; kwargs...) + +Keyword argument parser returning the appropriate `ProjectorAlgorithm` algorithm struct. +""" +function ProjectorAlgorithm(; + alg = Defaults.projector_alg, + decomposition_alg = (;), + trunc = (;), + verbosity = Defaults.projector_verbosity, + ) + # replace symbol with projector alg type + haskey(PROJECTOR_SYMBOLS, alg) || + throw(ArgumentError("unknown projector algorithm: $alg")) + alg_type = PROJECTOR_SYMBOLS[alg] + + # parse SVD forward & rrule algorithm + + decomposition_algorithm = if alg in [:HalfInfiniteProjector, :FullInfiniteProjector] + _alg_or_nt(SVDAdjoint, decomposition_alg) + elseif alg in [:C4vEighProjector] + _alg_or_nt(EighAdjoint, decomposition_alg) + elseif alg in [:C4vQRProjector] + _alg_or_nt(QRAdjoint, decomposition_alg) + end # TODO: how do we solve this in a proper way? + + # qr-ctmrg does not need truncation or degeneracy checks + if alg in [:C4vQRProjector] + return alg_type(decomposition_algorithm) + end + + # parse truncation scheme + truncation_strategy = if trunc isa TruncationStrategy + trunc + elseif trunc isa NamedTuple + _TruncationStrategy(; trunc...) + else + throw(ArgumentError("unknown trunc $trunc")) + end + + return alg_type(decomposition_algorithm, truncation_strategy, verbosity) +end + +""" + decomposition_algorithm(alg::ProjectorAlgorithm) + +Return the tensor decomposition algorithm of the `alg` projector algorithm. +""" +decomposition_algorithm(alg::ProjectorAlgorithm) = alg.decomposition_alg + +function truncation_strategy(alg::ProjectorAlgorithm, edge) + if alg.trunc isa FixedSpaceTruncation + tspace = space(edge, 1) + return isdual(tspace) ? truncspace(flip(tspace)) : truncspace(tspace) + else + return alg.trunc + end +end + +""" + _set_truncation(alg::ProjectorAlgorithm, trunc::TruncationStrategy) + +Update the truncation strategy of a given projector algorithm, keeping all other settings +the same. +""" +function _set_truncation(alg::ProjectorAlgorithm, trunc::TruncationStrategy) + alg = @set alg.trunc = trunc + return alg +end + +""" + _set_decomposition_truncation(alg::ProjectorAlgorithm, trunc::TruncationStrategy) + +Set the truncation strategy of the decomposition algorithm within a given projector +algorithm, keeping all other settings the same. +""" +function _set_decomposition_truncation(alg::ProjectorAlgorithm, trunc::TruncationStrategy) + decomp_alg = decomposition_algorithm(alg) + truncated_fwd_alg = TruncatedAlgorithm(decomp_alg.fwd_alg, trunc) + decomp_alg = @set decomp_alg.fwd_alg = truncated_fwd_alg + alg = @set alg.decomposition_alg = decomp_alg + return alg +end