Merge remote-tracking branch 'upstream/master' into bondenv-ntu

Author: Yue-Zhengyuan

src/algorithms/contractions/ctmrg/renormalize_corner.jl

@@ -336,7 +336,7 @@ Alternatively, provide the constituent tensors and perform the complete contract
 ```
 """
 function renormalize_southwest_corner((row, col), enlarged_env, P_left, P_right)
-    return renormalize_corner(
+    return renormalize_southwest_corner(
         enlarged_env[SOUTHWEST, row, col],
         P_left[WEST, row, col],
         P_right[SOUTH, row, _next(col, end)],
@@ -345,7 +345,7 @@ end
 function renormalize_southwest_corner(
         quadrant::AbstractTensorMap{T, S, N, N}, P_left, P_right
     ) where {T, S, N}
-    return renormalize_southwest_corner(quadrant, P_left, P_right)
+    return renormalize_corner(quadrant, P_left, P_right)
 end
 function renormalize_southwest_corner(
         E_south, C_southwest, E_west, P_left, P_right, A::PEPSSandwich

test/ctmrg/contractions.jl

@@ -0,0 +1,123 @@
+using Test
+using Random
+using PEPSKit
+using TensorKit
+
+using PEPSKit: eachcoordinate
+using PEPSKit: EnlargedCorner, simultaneous_projectors
+using PEPSKit: renormalize_northwest_corner, renormalize_northeast_corner, renormalize_southeast_corner, renormalize_southwest_corner
+
+# settings
+Random.seed!(91283219347)
+stype = ComplexF64
+ctm_alg = SimultaneousCTMRG(; projector_alg = :halfinfinite)
+
+function test_contractions(
+        Pspaces, Nspaces, Espaces, chis_north, chis_east, chis_south, chis_west,
+    )
+    peps = InfinitePEPS(randn, stype, Pspaces, Nspaces, Espaces)
+    env = CTMRGEnv(randn, stype, peps, chis_north, chis_east, chis_south, chis_west)
+
+    n = InfiniteSquareNetwork(peps)
+    coordinates = eachcoordinate(n)
+    dirs_and_coordinates = eachcoordinate(n, 1:4)
+
+    # initialize dense and sparse enlarged corners
+    sparse_enlarged_corners = map(dirs_and_coordinates) do co
+        return EnlargedCorner(n, env, co)
+    end
+    dense_enlarged_corners = map(TensorMap, sparse_enlarged_corners)
+
+    # compute projectors (doesn't matter how)
+    (P_left, P_right), info = simultaneous_projectors(
+        dense_enlarged_corners, env, ctm_alg.projector_alg
+    )
+
+    # test corner renormalization
+    return foreach(coordinates) do (r, c)
+        for renormalize_f in (
+                renormalize_northwest_corner, renormalize_northeast_corner,
+                renormalize_southeast_corner, renormalize_southwest_corner,
+            )
+            C_sparse = renormalize_f((r, c), sparse_enlarged_corners, P_left, P_right)
+            C_dense = renormalize_f((r, c), dense_enlarged_corners, P_left, P_right)
+            @test C_sparse ≈ C_dense
+        end
+    end
+
+    # TODO: test all other uncovered contracitions
+end
+
+@testset "Random Cartesian spaces" begin
+    unitcell = (3, 3)
+
+    Pspaces = ComplexSpace.(rand(2:3, unitcell...))
+    Nspaces = ComplexSpace.(rand(2:4, unitcell...))
+    Espaces = ComplexSpace.(rand(2:4, unitcell...))
+    chis_north = ComplexSpace.(rand(5:10, unitcell...))
+    chis_east = ComplexSpace.(rand(5:10, unitcell...))
+    chis_south = ComplexSpace.(rand(5:10, unitcell...))
+    chis_west = ComplexSpace.(rand(5:10, unitcell...))
+
+    test_contractions(
+        Pspaces, Nspaces, Espaces, chis_north, chis_east, chis_south, chis_west,
+    )
+end
+
+@testset "Specific U1 spaces" begin
+    unitcell = (2, 2)
+
+    PA = U1Space(-1 => 1, 0 => 1)
+    PB = U1Space(0 => 1, 1 => 1)
+    Vpeps = U1Space(-1 => 2, 0 => 1, 1 => 2)
+    Venv = U1Space(-2 => 2, -1 => 3, 0 => 4, 1 => 3, 2 => 2)
+
+    Pspaces = [PA PB; PB PA]
+    Nspaces = [Vpeps Vpeps'; Vpeps' Vpeps]
+    chis = [Venv Venv; Venv Venv]
+
+    test_contractions(Pspaces, Nspaces, Nspaces, chis, chis, chis, chis)
+
+    # 4x4 unit cell with all 32 inequivalent bonds
+    #
+    #    10     4     7    32
+    #     |     |     |     |
+    #  3--A--1--B--5--C--8--D--3
+    #     |     |     |     |
+    #     2     6     9    11
+    #     |     |     |     |
+    # 14--E-12--F-15--G-17--H-14
+    #     |     |     |     |
+    #    13    16    18    19
+    #     |     |     |     |
+    # 22--I-20--J-23--K-25--L-22
+    #     |     |     |     |
+    #    21    24    26    27
+    #     |     |     |     |
+    # 29--M-28--N-30--O-31--P-29
+    #     |     |     |     |
+    #    10     4     7    32
+
+    phys_space = Vect[U1Irrep](1 => 1, -1 => 1)
+    corner_space = Vect[U1Irrep](0 => 1, 1 => 1, -1 => 1)
+    vspaces = map(i -> Vect[U1Irrep](0 => 1 + i % 4, 1 => i ÷ 4 % 4, -2 => i ÷ 16), 1:32)
+    @test length(Set(vspaces)) == 32
+
+    Espaces = [
+        vspaces[1] vspaces[5] vspaces[8] vspaces[3]
+        vspaces[12] vspaces[15] vspaces[17] vspaces[14]
+        vspaces[20] vspaces[23] vspaces[25] vspaces[22]
+        vspaces[28] vspaces[30] vspaces[31] vspaces[29]
+    ]
+
+    Nspaces = [
+        vspaces[10] vspaces[4] vspaces[7] vspaces[32]
+        vspaces[2] vspaces[6] vspaces[9] vspaces[11]
+        vspaces[13] vspaces[16] vspaces[18] vspaces[19]
+        vspaces[21] vspaces[24] vspaces[26] vspaces[27]
+    ]
+    Pspaces = fill(phys_space, (4, 4))
+    chis = fill(corner_space, (4, 4))
+
+    test_contractions(Pspaces, Nspaces, Nspaces, chis, chis, chis, chis)
+end

test/examples/heisenberg.jl

@@ -55,7 +55,7 @@ end
 @testset "C4v AD optimization with scalartype T=$T and projector_alg=$projector_alg" for (T, projector_alg) in
     Iterators.product([Float64, ComplexF64], [:c4v_eigh, :c4v_qr])
     # initialize symmetric states
-    Random.seed!(1234567)
+    Random.seed!(123456789)
     symm = RotateReflect()
     H′ = heisenberg_XYZ_c4v(InfiniteSquare())
     H = T <: Real ? real(H′) : H′

test/runtests.jl

@@ -55,6 +55,9 @@ end
         @time @safetestset "correlation length" begin
             include("ctmrg/correlation_length.jl")
         end
+        @time @safetestset "Contractions" begin
+            include("ctmrg/contractions.jl")
+        end
     end
     if GROUP == "ALL" || GROUP == "BP"
         @time @safetestset "Unit cell bond matching" begin

test/timeevol/tf_ising_finiteT.jl

@@ -44,12 +44,13 @@ dt, nstep = 1.0e-3, 400
 β = dt * nstep
 
 # when g = 2, β = 0.4 and 2β = 0.8 belong to two phases (without and with nonzero σᶻ)
-@testset "Finite-T SU (bipartite = $(bipartite))" for bipartite in (true, false)
+@testset "Finite-T SU (force_mpo = $(force_mpo))" for force_mpo in (false, true)
     # use second order Trotter decomposition
     symmetrize_gates = true
+    bipartite = true
 
     # PEPO approach: results at β, or T = 2.5
-    alg = SimpleUpdate(; trunc = trunc_pepo, purified = false, bipartite)
+    alg = SimpleUpdate(; trunc = trunc_pepo, purified = false, bipartite, force_mpo)
     pepo, wts, info = time_evolve(pepo0, ham, dt, nstep, alg, wts0; symmetrize_gates)
 
     ## BP gauge fixing
@@ -72,7 +73,7 @@ dt, nstep = 1.0e-3, 400
     @test isapprox(abs.(result_2β), bm_2β, rtol = 1.0e-4)
 
     # Purification approach: results at 2β, or T = 1.25
-    alg = SimpleUpdate(; trunc = trunc_pepo, purified = true, bipartite)
+    alg = SimpleUpdate(; trunc = trunc_pepo, purified = true, bipartite, force_mpo)
     pepo, wts, info = time_evolve(pepo0, ham, dt, 2 * nstep, alg, wts0; symmetrize_gates)
     env = converge_env(InfinitePEPS(pepo), 8)
     result_2β′ = measure_mag(pepo, env; purified = true)