WIP

src/projmps.jl

@@ -214,10 +214,10 @@ function makesitediagonal(projmps::ProjMPS, sites::AbstractVector{Index})
 end
 
 function makesitediagonal(projmps::ProjMPS, tag::String)
-    mps_diagonal = Quantics.makesitediagonal(MPS(projmps), tag)
+    mps_diagonal = makesitediagonal(MPS(projmps), tag)
     projmps_diagonal = ProjMPS(mps_diagonal)
 
-    target_sites = Quantics.findallsiteinds_by_tag(
+    target_sites = findallsiteinds_by_tag(
         unique(ITensors.noprime.(Iterators.flatten(siteinds(projmps)))); tag=tag
     )
 
@@ -264,12 +264,9 @@ function extractdiagonal(
     return ProjMPS(MPS(tensors), projector)
 end
 
-function extractdiagonal(projmps::ProjMPS, site::Index{IndsT}) where {IndsT}
-    return Quantics.extractdiagonal(projmps, [site])
-end
 
 function extractdiagonal(projmps::ProjMPS, tag::String)::ProjMPS
-    targetsites = Quantics.findallsiteinds_by_tag(
+    targetsites = findallsiteinds_by_tag(
         unique(ITensors.noprime.(ProjMPSs._allsites(projmps))); tag=tag
     )
     return extractdiagonal(projmps, targetsites)

src/util.jl

@@ -116,3 +116,31 @@ function _extract_diagonal(t, site::Index{T}, site2::Index{T}) where {T<:Number}
     end
     return ITensor(newdata, restinds..., site)
 end
+
+
+"""
+Makes an MPS/MPO diagonal for a specified a site index `s`.
+On return, the data will be deep copied and the target core tensor will be diagonalized with an additional site index `s'`.
+"""
+function makesitediagonal(M::AbstractMPS, site::Index{T})::MPS where {T}
+    M_ = deepcopy(MPO(collect(M)))
+
+    target_site::Int = only(findsites(M_, site))
+    M_[target_site] = _asdiagonal(M_[target_site], site)
+
+    return MPS(collect(M_))
+end
+
+function makesitediagonal(M::AbstractMPS, tag::String)::MPS
+    M_ = deepcopy(MPO(collect(M)))
+    sites = siteinds(M_)
+
+    target_positions = findallsites_by_tag(siteinds(M_); tag=tag)
+
+    for t in eachindex(target_positions)
+        i, j = target_positions[t]
+        M_[i] = _asdiagonal(M_[i], sites[i][j])
+    end
+
+    return MPS(collect(M_))
+end

test/runtests.jl

@@ -6,10 +6,11 @@ using Test
 
 include("_util.jl")
 
-include("projector_tests.jl")
-include("projmps_tests.jl")
-include("blockedmps_tests.jl")
-include("contract_tests.jl")
-include("patching_tests.jl")
+#include("projector_tests.jl")
+#include("projmps_tests.jl")
+#include("blockedmps_tests.jl")
+#include("contract_tests.jl")
+#include("patching_tests.jl")
+include("util_tests.jl")
 
-include("automul_tests.jl")
+#include("automul_tests.jl")

test/util_tests.jl

@@ -0,0 +1,81 @@
+using Test
+
+using ITensors
+import ProjMPSs: ProjMPSs, Projector, project, ProjMPS, projcontract, BlockedMPS, rearrange_siteinds, makesitediagonal, extractdiagonal
+import FastMPOContractions as FMPOC
+
+@testset "util.jl" begin
+    @testset "rearrange_siteinds" begin
+        N = 3
+        sitesx = [Index(2, "x=$n") for n in 1:N]
+        sitesy = [Index(2, "y=$n") for n in 1:N]
+        sitesz = [Index(2, "z=$n") for n in 1:N]
+        sites = collect(collect.(zip(sitesx, sitesy, sitesz)))
+
+        Ψ = MPS(collect(_random_mpo(sites)))
+
+        prjΨ = ProjMPS(Ψ)
+        prjΨ1 = project(prjΨ, Dict(sitesx[1] => 1))
+
+        sitesxy = collect(collect.(zip(sitesx, sitesy)))
+        sites_rearranged = Vector{Index{Int}}[]
+        for i in 1:N
+            push!(sites_rearranged, sitesxy[i])
+            push!(sites_rearranged, [sitesz[i]])
+        end
+        prjΨ1_rearranged = rearrange_siteinds(prjΨ1, sites_rearranged)
+
+        @test reduce(*, MPS(prjΨ1)) ≈ reduce(*, MPS(prjΨ1_rearranged))
+        @test ProjMPSs.siteinds(prjΨ1_rearranged) == sites_rearranged
+    end
+
+    @testset "makesitediagonal and extractdiagonal" begin
+        N = 3
+        sitesx = [Index(2, "x=$n") for n in 1:N]
+        sitesy = [Index(2, "y=$n") for n in 1:N]
+        sitesz = [Index(2, "z=$n") for n in 1:N]
+
+        sitesxy_vec = [[x, y] for (x, y) in zip(sitesx, sitesy)]
+        sitesz_vec = [[z] for z in sitesz]
+        sites = [x for pair in zip(sitesxy_vec, sitesz_vec) for x in pair]
+
+        Ψ = MPS(collect(_random_mpo(sites)))
+
+        prjΨ = ProjMPS(Ψ)
+        prjΨ1 = project(prjΨ, Dict(sitesx[1] => 1))
+
+        prjΨ1_diagonalz = makesitediagonal(prjΨ1, "y")
+        sites_diagonalz = Iterators.flatten(siteinds(prjΨ1_diagonalz))
+
+        psi_diag = prod(prjΨ1_diagonalz.data)
+        psi = prod(prjΨ1.data)
+
+        @test extractdiagonal(prjΨ1_diagonalz, "y") ≈ prjΨ1
+
+        for indval in eachindval(sites_diagonalz...)
+            ind = first.(indval)
+            val = last.(indval)
+
+            index_dict = Dict{Index{Int},Vector{Int}}()
+            for (i, el) in enumerate(ind)
+                baseind = noprime(el)
+                if haskey(index_dict, baseind)
+                    push!(index_dict[baseind], i)
+                else
+                    index_dict[baseind] = [i]
+                end
+            end
+            repeated_indices = [is for is in values(index_dict) if length(is) > 1]
+
+            isdiagonalelement = all(allequal(val[i] for i in is)
+            for is in repeated_indices)
+
+            if isdiagonalelement
+                nondiaginds = unique(noprime(i) => v for (i, v) in indval)
+                @test psi_diag[indval...] == psi[nondiaginds...]
+            else
+                @test iszero(psi_diag[indval...])
+            end
+        end
+    end
+end