Skip to content

update left and right_virtualspace #205

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
lkdvos merged 1 commit into from
Dec 12, 2024
Merged

update left and right_virtualspace #205

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
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: 0 additions & 2 deletions src/MPSKit.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -60,8 +60,6 @@ export exact_diagonalization
export TransferMatrix
export transfer_left, transfer_right

@deprecate virtualspace left_virtualspace # there is a possible ambiguity when C isn't square, necessitating specifying left or right virtualspace

# Abstract type defs
abstract type Algorithm end

Expand Down
6 changes: 3 additions & 3 deletions src/algorithms/expval.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -66,8 +66,8 @@ function expectation_value(ψ::AbstractMPS, (inds, O)::Pair)
# left side
T = storagetype(site_type(ψ))
@plansor Vl[-1 -2; -3] := isomorphism(T,
left_virtualspace(ψ, sites[1] - 1),
left_virtualspace(ψ, sites[1] - 1))[-1; -3] *
left_virtualspace(ψ, sites[1]),
left_virtualspace(ψ, sites[1]))[-1; -3] *
conj(Ut[-2])

# middle
Expand Down Expand Up @@ -104,7 +104,7 @@ function expectation_value(ψ::InfiniteMPS, H::InfiniteMPOHamiltonian,
envs::AbstractMPSEnvironments=environments(ψ, H))
return sum(1:length(ψ)) do i
util = fill_data!(similar(ψ.AL[1], right_virtualspace(H, i)[end]), one)
@plansor GR[-1 -2; -3] := r_LL(ψ, i)[-1; -3] * conj(util[-2])
@plansor GR[-1 -2; -3] := r_LL(ψ, i)[-1; -3] * util[-2]
return contract_mpo_expval(ψ.AL[i], leftenv(envs, i, ψ), H[i][:, 1, 1, end], GR)
end
end
Expand Down
2 changes: 1 addition & 1 deletion src/algorithms/fidelity_susceptibility.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@
Tos = LeftGaugedQP(rand, state)
for (i, ac) in enumerate(state.AC)
temp = ∂∂AC(i, state, V, venvs) * ac
help = fill_data!(similar(ac, utilleg(Tos)), one)
help = fill_data!(similar(ac, auxiliaryspace(Tos)), one)

Check warning on line 14 in src/algorithms/fidelity_susceptibility.jl

View check run for this annotation

Codecov / codecov/patch

src/algorithms/fidelity_susceptibility.jl#L14 

Added line #L14 was not covered by tests
@plansor Tos[i][-1 -2; -3 -4] := temp[-1 -2; -4] * help[-3]
end

Expand Down
2 changes: 1 addition & 1 deletion src/algorithms/timestep/timeevmpo.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -327,7 +327,7 @@ function make_time_mpo(H::InfiniteMPOHamiltonian{T}, dt, alg::WII) where {T}
Vᵣ = right_virtualspace(H, i)[1:(end - 1)]
P = physicalspace(H, i)

h′ = similar(H[i], Vₗ ⊗ P ← P ⊗ Vᵣ')
h′ = similar(H[i], Vₗ ⊗ P ← P ⊗ Vᵣ)
h′[2:end, 1, 1, 2:end] = WA[i]
h′[2:end, 1, 1, 1] = WB[i]
h′[1, 1, 1, 2:end] = WC[i]
Expand Down
22 changes: 11 additions & 11 deletions src/algorithms/toolbox.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -63,10 +63,10 @@
"""
function transfer_spectrum(above::InfiniteMPS; below=above, tol=Defaults.tol, num_vals=20,
sector=first(sectors(oneunit(left_virtualspace(above, 1)))))
init = randomize!(similar(above.AL[1], left_virtualspace(below, 0),
ℂ[typeof(sector)](sector => 1)' * left_virtualspace(above, 0)))
init = randomize!(similar(above.AL[1], left_virtualspace(below, 1),
ℂ[typeof(sector)](sector => 1)' * left_virtualspace(above, 1)))

transferspace = fuse(left_virtualspace(above, 0) * left_virtualspace(below, 0)')
transferspace = fuse(left_virtualspace(above, 1) * left_virtualspace(below, 1)')
num_vals = min(dim(transferspace, sector), num_vals) # we can ask at most this many values
eigenvals, eigenvecs, convhist = eigsolve(flip(TransferMatrix(above.AL, below.AL)),
init, num_vals, :LM; tol=tol)
Expand Down Expand Up @@ -239,26 +239,26 @@

# allocate output
output = Vector{O}(undef, L)
V_wrap = left_virtualspace(mpo, 1)'
V_wrap = left_virtualspace(mpo, 1)

Check warning on line 242 in src/algorithms/toolbox.jl

View check run for this annotation

Codecov / codecov/patch

src/algorithms/toolbox.jl#L242 

Added line #L242 was not covered by tests
ST = storagetype(O)

util = fill!(similar(mpo[1], oneunit(V_wrap)), one(scalartype(O)))
@plansor cup[-1; -2 -3] := id(ST, V_wrap)[-3; -2] * util[-1]
@plansor cup[-1; -2 -3] := id(ST, V_wrap)[-2; -3] * util[-1]

Check warning on line 246 in src/algorithms/toolbox.jl

View check run for this annotation

Codecov / codecov/patch

src/algorithms/toolbox.jl#L246 

Added line #L246 was not covered by tests

local F_right
for i in 1:L
V_left = i == 1 ? oneunit(V_wrap) : fuse(V_wrap ⊗ left_virtualspace(mpo, i))
V_right = i == L ? oneunit(V_wrap) : fuse(V_wrap' ⊗ right_virtualspace(mpo, i)')
V_right = i == L ? oneunit(V_wrap) : fuse(V_wrap ⊗ right_virtualspace(mpo, i))

Check warning on line 251 in src/algorithms/toolbox.jl

View check run for this annotation

Codecov / codecov/patch

src/algorithms/toolbox.jl#L251 

Added line #L251 was not covered by tests
output[i] = similar(mpo[i],
V_left * physicalspace(mpo, i) ←
physicalspace(mpo, i) * V_right)
F_left = i == 1 ? cup : F_right
F_right = i == L ? cup :
isomorphism(ST, V_right ← V_wrap * right_virtualspace(mpo, i)')
@plansor output[i][-1 -2; -3 -4] = F_left[-1; 1 2] *
τ[-3 1; 4 3] *
mpo[i][2 -2; 3 5] *
conj(F_right[-4; 4 5])
isomorphism(ST, V_right ← V_wrap' * right_virtualspace(mpo, i))
@plansor contractcheck = true output[i][-1 -2; -3 -4] = F_left[-1; 1 2] *

Check warning on line 258 in src/algorithms/toolbox.jl

View check run for this annotation

Codecov / codecov/patch

src/algorithms/toolbox.jl#L258 

Added line #L258 was not covered by tests
τ[-3 1; 4 3] *
mpo[i][2 -2; 3 5] *
conj(F_right[-4; 4 5])
end

mpo isa SparseMPO && dropzeros!.(output) # the above process fills sparse mpos with zeros.
Expand Down
16 changes: 7 additions & 9 deletions src/environments/abstract_envs.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -14,12 +14,10 @@ Base.unlock(envs::AbstractMPSEnvironments) = unlock(envs.lock);
# Allocating tensors
# ------------------

# TODO: fix the fucking left/right virtualspace bullshit
# TODO: storagetype stuff
function allocate_GL(bra::AbstractMPS, mpo::AbstractMPO, ket::AbstractMPS, i::Int)
T = Base.promote_type(scalartype(bra), scalartype(mpo), scalartype(ket))
V = left_virtualspace(bra, i - 1) ⊗ left_virtualspace(mpo, i)' ←
left_virtualspace(ket, i - 1)
V = left_virtualspace(bra, i) ⊗ left_virtualspace(mpo, i)' ←
left_virtualspace(ket, i)
if V isa BlockTensorKit.TensorMapSumSpace
TT = blocktensormaptype(spacetype(bra), numout(V), numin(V), T)
else
Expand All @@ -30,7 +28,7 @@ end

function allocate_GR(bra::AbstractMPS, mpo::AbstractMPO, ket::AbstractMPS, i::Int)
T = Base.promote_type(scalartype(bra), scalartype(mpo), scalartype(ket))
V = right_virtualspace(ket, i) ⊗ right_virtualspace(mpo, i)' ←
V = right_virtualspace(ket, i) ⊗ right_virtualspace(mpo, i) ←
right_virtualspace(bra, i)
if V isa BlockTensorKit.TensorMapSumSpace
TT = blocktensormaptype(spacetype(bra), numout(V), numin(V), T)
Expand All @@ -42,8 +40,8 @@ end

function allocate_GBL(bra::QP, mpo::AbstractMPO, ket::QP, i::Int)
T = Base.promote_type(scalartype(bra), scalartype(mpo), scalartype(ket))
V = left_virtualspace(bra.left_gs, i - 1) ⊗ left_virtualspace(mpo, i)' ←
auxiliaryspace(ket)' ⊗ left_virtualspace(ket.left_gs, i - 1)
V = left_virtualspace(bra, i) ⊗ left_virtualspace(mpo, i)' ←
auxiliaryspace(ket)' ⊗ left_virtualspace(ket, i)
if V isa BlockTensorKit.TensorMapSumSpace
TT = blocktensormaptype(spacetype(bra), numout(V), numin(V), T)
else
Expand All @@ -54,8 +52,8 @@ end

function allocate_GBR(bra::QP, mpo::AbstractMPO, ket::QP, i::Int)
T = Base.promote_type(scalartype(bra), scalartype(mpo), scalartype(ket))
V = right_virtualspace(ket.right_gs, i) ⊗ right_virtualspace(mpo, i)' ←
auxiliaryspace(ket)' ⊗ right_virtualspace(bra.right_gs, i)
V = right_virtualspace(ket, i) ⊗ right_virtualspace(mpo, i) ←
auxiliaryspace(ket)' ⊗ right_virtualspace(bra, i)
if V isa BlockTensorKit.TensorMapSumSpace
TT = blocktensormaptype(spacetype(bra), numout(V), numin(V), T)
else
Expand Down
18 changes: 3 additions & 15 deletions src/environments/finite_envs.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -34,30 +34,18 @@ function environments(below, operator, above, leftstart, rightstart)
rightenvs)
end

function environments(below::FiniteMPS{S}, O::DenseMPO, above=nothing) where {S}
N = length(below)
leftstart = isomorphism(storagetype(S),
left_virtualspace(below, 0) ⊗ space(O[1], 1)' ←
left_virtualspace(something(above, below), 0))
rightstart = isomorphism(storagetype(S),
right_virtualspace(something(above, below), N) ⊗
space(O[N], 4)' ←
right_virtualspace(below, length(below)))
return environments(below, O, above, leftstart, rightstart)
end

function environments(below::FiniteMPS{S}, O::Union{FiniteMPO,FiniteMPOHamiltonian},
above=nothing) where {S}
Vl_bot = left_virtualspace(below, 0)
Vl_bot = left_virtualspace(below, 1)
Vl_mid = left_virtualspace(O, 1)
Vl_top = isnothing(above) ? left_virtualspace(below, 0) : left_virtualspace(above, 0)
Vl_top = isnothing(above) ? left_virtualspace(below, 1) : left_virtualspace(above, 1)
leftstart = isomorphism(storagetype(S), Vl_bot ⊗ Vl_mid' ← Vl_top)

N = length(below)
Vr_bot = right_virtualspace(below, N)
Vr_mid = right_virtualspace(O, N)
Vr_top = isnothing(above) ? right_virtualspace(below, N) : right_virtualspace(above, N)
rightstart = isomorphism(storagetype(S), Vr_top ⊗ Vr_mid' ← Vr_bot)
rightstart = isomorphism(storagetype(S), Vr_top ⊗ Vr_mid ← Vr_bot)

return environments(below, O, above, leftstart, rightstart)
end
Expand Down
6 changes: 3 additions & 3 deletions src/operators/abstractmpo.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -156,7 +156,7 @@ Compute the mpo tensor that arises from multiplying MPOs.
function fuse_mul_mpo(O1::MPOTensor, O2::MPOTensor)
T = promote_type(scalartype(O1), scalartype(O2))
F_left = fuser(T, left_virtualspace(O2), left_virtualspace(O1))
F_right = fuser(T, right_virtualspace(O2)', right_virtualspace(O1)')
F_right = fuser(T, right_virtualspace(O2), right_virtualspace(O1))
@plansor O[-1 -2; -3 -4] := F_left[-1; 1 2] *
O2[1 5; -3 3] *
O1[2 -2; 5 4] *
Expand Down Expand Up @@ -203,7 +203,7 @@ function add_physical_charge(O::BraidingTensor, charge::Sector)
sectortype(O) === typeof(charge) || throw(SectorMismatch())
auxspace = Vect[typeof(charge)](charge => 1)
V = left_virtualspace(O) ⊗ fuse(physicalspace(O), auxspace) ←
fuse(physicalspace(O), auxspace) ⊗ right_virtualspace(O)'
fuse(physicalspace(O), auxspace) ⊗ right_virtualspace(O)
return BraidingTensor{scalartype(O)}(V)
end
function add_physical_charge(O::AbstractBlockTensorMap{<:Any,<:Any,2,2}, charge::Sector)
Expand All @@ -213,7 +213,7 @@ function add_physical_charge(O::AbstractBlockTensorMap{<:Any,<:Any,2,2}, charge:

Odst = similar(O,
left_virtualspace(O) ⊗ fuse(physicalspace(O), auxspace) ←
fuse(physicalspace(O), auxspace) ⊗ right_virtualspace(O)')
fuse(physicalspace(O), auxspace) ⊗ right_virtualspace(O))
for (I, v) in nonzero_pairs(O)
Odst[I] = add_physical_charge(v, charge)
end
Expand Down
55 changes: 23 additions & 32 deletions src/operators/mpo.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,7 @@ const FiniteMPO{O<:MPOTensor} = MPO{O,Vector{O}}

function FiniteMPO(Os::AbstractVector{O}) where {O<:MPOTensor}
for i in eachindex(Os)[1:(end - 1)]
dual(right_virtualspace(Os[i])) == left_virtualspace(Os[i + 1]) ||
right_virtualspace(Os[i]) == left_virtualspace(Os[i + 1]) ||
throw(SpaceMismatch("unmatching virtual spaces at site $i"))
end
return FiniteMPO{O}(Os)
Expand All @@ -38,7 +38,7 @@ const InfiniteMPO{O<:MPOTensor} = MPO{O,PeriodicVector{O}}

function InfiniteMPO(Os::AbstractVector{O}) where {O<:MPOTensor}
for i in eachindex(Os)
dual(right_virtualspace(Os[i])) == left_virtualspace(Os[mod1(i + 1, end)]) ||
right_virtualspace(Os[i]) == left_virtualspace(Os[mod1(i + 1, end)]) ||
throw(SpaceMismatch("umatching virtual spaces at site $i"))
end
return InfiniteMPO{O}(Os)
Expand Down Expand Up @@ -139,8 +139,8 @@ function Base.convert(::Type{TensorMap}, mpo::FiniteMPO)
U_left = ones(scalartype(mpo), V_left)'

V_right = right_virtualspace(mpo, length(mpo))
@assert V_right == oneunit(V_right)'
U_right = ones(scalartype(mpo), V_right')
@assert V_right == oneunit(V_right)
U_right = ones(scalartype(mpo), V_right)

tensors = vcat(U_left, parent(mpo), U_right)
indices = [[i, -i, -(2N - i + 1), i + 1] for i in 1:length(mpo)]
Expand All @@ -166,10 +166,10 @@ function Base.:+(mpo1::FiniteMPO{TO}, mpo2::FiniteMPO{TO}) where {TO}
A = storagetype(TO)

# left half
F₁ = isometry(A, (right_virtualspace(mpo1, 1) ⊕ right_virtualspace(mpo2, 1))',
right_virtualspace(mpo1, 1)')
F₁ = isometry(A, (right_virtualspace(mpo1, 1) ⊕ right_virtualspace(mpo2, 1)),
right_virtualspace(mpo1, 1))
F₂ = leftnull(F₁)
@assert _lastspace(F₂) == right_virtualspace(mpo2, 1)
@assert _lastspace(F₂) == right_virtualspace(mpo2, 1)'

@plansor O[-3 -1 -2; -4] := mpo1[1][-1 -2; -3 1] * conj(F₁[-4; 1]) +
mpo2[1][-1 -2; -3 1] * conj(F₂[-4; 1])
Expand All @@ -184,10 +184,10 @@ function Base.:+(mpo1::FiniteMPO{TO}, mpo2::FiniteMPO{TO}) where {TO}
@plansor O₂[-1 -2; -3 -4] := R[-1; 1] * F₂[1; 2] * mpo2[i][2 -2; -3 -4]

# incorporate fusers from right side
F₁ = isometry(A, (right_virtualspace(mpo1, i) ⊕ right_virtualspace(mpo2, i))',
right_virtualspace(mpo1, i)')
F₁ = isometry(A, (right_virtualspace(mpo1, i) ⊕ right_virtualspace(mpo2, i)),
right_virtualspace(mpo1, i))
F₂ = leftnull(F₁)
@assert _lastspace(F₂) == right_virtualspace(mpo2, i)
@assert _lastspace(F₂) == right_virtualspace(mpo2, i)'
@plansor O[-3 -1 -2; -4] := O₁[-1 -2; -3 1] * conj(F₁[-4; 1]) +
O₂[-1 -2; -3 1] * conj(F₂[-4; 1])

Expand Down Expand Up @@ -248,8 +248,8 @@ function Base.:*(mpo1::FiniteMPO{TO}, mpo2::FiniteMPO{TO}) where {TO}
S = spacetype(TO)
if (left_virtualspace(mpo1, 1) != oneunit(S) ||
left_virtualspace(mpo2, 1) != oneunit(S)) ||
(right_virtualspace(mpo1, N)' != oneunit(S) ||
right_virtualspace(mpo2, N)' != oneunit(S))
(right_virtualspace(mpo1, N) != oneunit(S) ||
right_virtualspace(mpo2, N) != oneunit(S))
@warn "left/right virtual space is not trivial, fusion may not be unique"
# this is a warning because technically any isomorphism that fuses the left/right
# would work and for now I dont feel like figuring out if this is important
Expand All @@ -261,13 +261,8 @@ function Base.:*(mpo1::FiniteMPO{TO}, mpo2::FiniteMPO{TO}) where {TO}
# note order of mpos: mpo1 * mpo2 * state -> mpo2 on top of mpo1
local Fᵣ # trick to make Fᵣ defined in the loop
for i in 1:N
Fₗ = i != 1 ? Fᵣ :
isomorphism(A, fuse(left_virtualspace(mpo2, i), left_virtualspace(mpo1, i)),
left_virtualspace(mpo2, i) * left_virtualspace(mpo1, i))
Fᵣ = isomorphism(A,
fuse(right_virtualspace(mpo2, i)', right_virtualspace(mpo1, i)'),
right_virtualspace(mpo2, i)' * right_virtualspace(mpo1, i)')

Fₗ = i != 1 ? Fᵣ : fuser(A, left_virtualspace(mpo2, i), left_virtualspace(mpo1, i))
Fᵣ = fuser(A, right_virtualspace(mpo2, i), right_virtualspace(mpo1, i))
@plansor O[i][-1 -2; -3 -4] := Fₗ[-1; 1 4] * mpo2[i][1 2; -3 3] *
mpo1[i][4 -2; 2 5] *
conj(Fᵣ[-4; 3 5])
Expand All @@ -284,11 +279,8 @@ function Base.:*(mpo::FiniteMPO, mps::FiniteMPS)

local Fᵣ # trick to make Fᵣ defined in the loop
for i in 1:length(mps)
Fₗ = i != 1 ? Fᵣ :
isomorphism(TT, fuse(left_virtualspace(A[i]), left_virtualspace(mpo, i)),
left_virtualspace(A[i]) * left_virtualspace(mpo, i))
Fᵣ = isomorphism(TT, fuse(right_virtualspace(A[i])', right_virtualspace(mpo, i)'),
right_virtualspace(A[i])' * right_virtualspace(mpo, i)')
Fₗ = i != 1 ? Fᵣ : fuser(TT, left_virtualspace(mps, i), left_virtualspace(mpo, i))
Fᵣ = fuser(TT, right_virtualspace(mps, i), right_virtualspace(mpo, i))
A[i] = _fuse_mpo_mps(mpo[i], A[i], Fₗ, Fᵣ)
end

Expand All @@ -298,12 +290,11 @@ function Base.:*(mpo::FiniteMPO, mps::FiniteMPS)
end

function Base.:*(mpo::InfiniteMPO, mps::InfiniteMPS)
check_length(mpo, mps)
L = check_length(mpo, mps)
T = promote_type(scalartype(mpo), scalartype(mps))
fusers = PeriodicArray(map(mps.AL, mpo) do al, mp
return fuser(T, _firstspace(al), _firstspace(mp))
end)
As = map(1:length(mps)) do i
fusers = PeriodicArray(fuser.(T, left_virtualspace.(Ref(mps), 1:L),
left_virtualspace.(Ref(mpo), 1:L)))
As = map(1:L) do i
return _fuse_mpo_mps(mpo[i], mps.AL[i], fusers[i], fusers[i + 1])
end
return changebonds(InfiniteMPS(As), SvdCut(; trscheme=notrunc()))
Expand Down Expand Up @@ -346,14 +337,14 @@ function TensorKit.dot(bra::FiniteMPS{T}, mpo::FiniteMPO, ket::FiniteMPS{T}) whe
Nhalf = N ÷ 2
# left half
ρ_left = isomorphism(storagetype(T),
left_virtualspace(bra, 0) ⊗ left_virtualspace(mpo, 1)',
left_virtualspace(ket, 0))
left_virtualspace(bra, 1) ⊗ left_virtualspace(mpo, 1)',
left_virtualspace(ket, 1))
T_left = TransferMatrix(ket.AL[1:Nhalf], mpo[1:Nhalf], bra.AL[1:Nhalf])
ρ_left = ρ_left * T_left

# right half
ρ_right = isomorphism(storagetype(T),
right_virtualspace(ket, N) ⊗ right_virtualspace(mpo, N)',
right_virtualspace(ket, N) ⊗ right_virtualspace(mpo, N),
right_virtualspace(ket, length(ket)))
T_right = TransferMatrix(ket.AR[(Nhalf + 1):end], mpo[(Nhalf + 1):end],
bra.AR[(Nhalf + 1):end])
Expand Down
Loading
Loading