Add iPEPS SUWeight to CTMRGEnv conversion

src/environments/suweight.jl

@@ -370,3 +370,49 @@ function Base.rot180(wts::SUWeight)
     wts_y = circshift(rot180(wts[2, :, :]), (1, 0))
     return SUWeight(wts_x, wts_y)
 end
+
+function _CTMRGEnv(wts::SUWeight, flips::Array{Bool, 3})
+    @assert size(wts) == size(flips)
+    _, Nr, Nc = size(wts)
+    edges = map(Iterators.product(1:4, 1:Nr, 1:Nc)) do (d, r, c)
+        wt_idx = if d == NORTH
+            CartesianIndex(2, _next(r, Nr), c)
+        elseif d == EAST
+            CartesianIndex(1, r, _prev(c, Nc))
+        elseif d == SOUTH
+            CartesianIndex(2, r, c)
+        else # WEST
+            CartesianIndex(1, r, c)
+        end
+        wt = deepcopy(wts[wt_idx])
+        if (flips[wt_idx] && d in (NORTH, EAST)) || (!flips[wt_idx] && d in (SOUTH, WEST))
+            wt = permute(wt, ((2,), (1,)))
+        end
+        twistdual!(wt, 1)
+        wt = insertrightunit(insertleftunit(wt, 1))
+        return permute(wt, ((1, 2, 3), (4,)))
+    end
+    corners = map(CartesianIndices(edges)) do idx
+        return TensorKit.id(Float64, codomain(edges[idx], 1))
+    end
+    return CTMRGEnv(corners, edges)
+end
+
+"""
+    CTMRGEnv(wts::SUWeight, peps::InfinitePEPS)
+
+Construct a CTMRG environment with bond dimension χ = 1 
+for an InfinitePEPS `peps` from the accompanying SUWeight `wts`.
+The scalartype of the returned environment is always `Float64`.
+"""
+function CTMRGEnv(wts::SUWeight, peps::InfinitePEPS)
+    _, Nr, Nc = size(wts)
+    flips = map(Iterators.product(1:2, 1:Nr, 1:Nc)) do (d, r, c)
+        return if d == 1
+            isdual(domain(peps[r, c], EAST))
+        else
+            isdual(domain(peps[r, c], NORTH))
+        end
+    end
+    return _CTMRGEnv(wts, flips)
+end

test/ctmrg/suweight.jl

@@ -0,0 +1,40 @@
+using Test
+using Random
+using TensorKit
+using PEPSKit
+using PEPSKit: str
+
+function rand_wts(peps::InfinitePEPS)
+    wts = SUWeight(peps)
+    for idx in CartesianIndices(wts.data)
+        n = length(wts.data[idx].data)
+        wts.data[idx].data[:] = sort(rand(Float64, n); lt = !isless)
+    end
+    return wts
+end
+
+function su_rdm_1x1(row::Int, col::Int, peps::InfinitePEPS, wts::SUWeight)
+    Nr, Nc = size(peps)
+    @assert 1 <= row <= Nr && 1 <= col <= Nc
+    t = peps.A[row, col]
+    t = absorb_weight(t, wts, row, col, Tuple(1:4))
+    ρ = t * t'
+    return ρ / str(ρ)
+end
+
+Vps = [Vect[U1Irrep](0 => 2, 1 => 2, -1 => 1), Vect[FermionParity](0 => 1, 1 => 2)]
+Vvs = [Vect[U1Irrep](0 => 3, 1 => 1, -1 => 2), Vect[FermionParity](0 => 2, 1 => 2)]
+
+for (Vp, Vv) in zip(Vps, Vvs)
+    Nspaces = [Vv Vv' Vv; Vv' Vv Vv']
+    Espaces = [Vv Vv Vv'; Vv Vv' Vv']
+    Pspaces = fill(Vp, size(Nspaces))
+    peps = InfinitePEPS(randn, ComplexF64, Pspaces, Nspaces, Espaces)
+    wts = rand_wts(peps)
+    env = CTMRGEnv(wts, peps)
+    for (r, c) in Tuple.(CartesianIndices(peps.A))
+        ρ1 = su_rdm_1x1(r, c, peps, wts)
+        ρ2 = reduced_densitymatrix(((r, c),), peps, env)
+        @test ρ1 ≈ ρ2
+    end
+end

test/runtests.jl

@@ -20,6 +20,9 @@ end
         @time @safetestset ":fixed CTMRG iteration scheme" begin
             include("ctmrg/fixed_iterscheme.jl")
         end
+        @time @safetestset "SUWeight conversion" begin
+            include("ctmrg/suweight.jl")
+        end
         @time @safetestset "Flavors" begin
             include("ctmrg/flavors.jl")
         end