Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
26 changes: 13 additions & 13 deletions src/algorithms/toolbox.jl
Original file line number Diff line number Diff line change
Expand Up @@ -36,8 +36,8 @@ function infinite_temperature_density_matrix(H::MPOHamiltonian)
end

"""
calc_galerkin(above, operator, below, envs)
calc_galerkin(pos, above, operator, below, envs)
calc_galerkin(below, operator, above, envs)
calc_galerkin(pos, below, operator, above, envs)

Calculate the Galerkin error, which is the error between the solution of the original problem, and the solution of the problem projected on the tangent space.
Concretely, this is the overlap of the current state with the single-site derivative, projected onto the nullspace of the current state:
Expand All @@ -46,25 +46,25 @@ Concretely, this is the overlap of the current state with the single-site deriva
\\epsilon = |VL * (VL' * \\frac{above}{\\partial AC_{pos}})|
```
"""
function calc_galerkin(pos::Int, above::Union{InfiniteMPS,FiniteMPS,WindowMPS}, operator,
below, envs)
AC´ = ∂∂AC(pos, above, operator, envs) * above.AC[pos]
function calc_galerkin(pos::Int, below::Union{InfiniteMPS,FiniteMPS,WindowMPS}, operator,
above, envs)
AC´ = ∂∂AC(pos, below, operator, envs) * above.AC[pos]
Comment thread
lkdvos marked this conversation as resolved.
normalize!(AC´)
out = add!(AC´, below.AL[pos] * below.AL[pos]' * AC´, -1)
return norm(out)
end
function calc_galerkin(pos::CartesianIndex{2}, above::MultilineMPS, operator::MultilineMPO,
below::MultilineMPS, envs::MultilineEnvironments)
function calc_galerkin(pos::CartesianIndex{2}, below::MultilineMPS, operator::MultilineMPO,
above::MultilineMPS, envs::MultilineEnvironments)
row, col = pos.I
return calc_galerkin(col, above[row], operator[row], below[row + 1], envs[row])
return calc_galerkin(col, below[row + 1], operator[row], above[row], envs[row])
end
function calc_galerkin(above::Union{InfiniteMPS,FiniteMPS,WindowMPS}, operator,
below, envs)
return maximum(pos -> calc_galerkin(pos, above, operator, below, envs), 1:length(above))
function calc_galerkin(below::Union{InfiniteMPS,FiniteMPS,WindowMPS}, operator,
above, envs)
return maximum(pos -> calc_galerkin(pos, below, operator, above, envs), 1:length(above))
end
function calc_galerkin(above::MultilineMPS, operator::MultilineMPO, below::MultilineMPS,
function calc_galerkin(below::MultilineMPS, operator::MultilineMPO, above::MultilineMPS,
envs::MultilineEnvironments)
return maximum(pos -> calc_galerkin(pos, above, operator, below, envs),
return maximum(pos -> calc_galerkin(pos, below, operator, above, envs),
CartesianIndices(size(above)))
end

Expand Down
4 changes: 2 additions & 2 deletions src/environments/abstract_envs.jl
Original file line number Diff line number Diff line change
Expand Up @@ -64,7 +64,7 @@ end
# Environment algorithms
# ----------------------
"""
environment_alg(above, operator, below; kwargs...)
environment_alg(below, operator, above; kwargs...)

Determine an appropriate algorithm for computing the environments, based on the given `kwargs...`.
"""
Expand All @@ -76,7 +76,7 @@ function environment_alg(::Union{InfiniteMPS,MultilineMPS},
eager=true)
return Arnoldi(; tol, maxiter, krylovdim, verbosity, eager)
end
function environment_alg(above, ::InfiniteMPOHamiltonian, below;
function environment_alg(below, ::InfiniteMPOHamiltonian, above;
tol=Defaults.tol, maxiter=Defaults.maxiter,
krylovdim=Defaults.krylovdim, verbosity=Defaults.VERBOSE_NONE)
return GMRES(; tol, maxiter, krylovdim, verbosity)
Expand Down
94 changes: 47 additions & 47 deletions src/environments/infinite_envs.jl
Original file line number Diff line number Diff line change
Expand Up @@ -18,19 +18,19 @@ Base.length(envs::InfiniteEnvironments) = length(envs.GLs)
leftenv(envs::InfiniteEnvironments, site::Int, state) = envs.GLs[site]
rightenv(envs::InfiniteEnvironments, site::Int, state) = envs.GRs[site]

function environments(above::InfiniteMPS,
function environments(below::InfiniteMPS,
operator::Union{InfiniteMPO,InfiniteMPOHamiltonian},
below::InfiniteMPS=above;
above::InfiniteMPS=below;
kwargs...)
GLs, GRs = initialize_environments(above, operator, below)
GLs, GRs = initialize_environments(below, operator, above)
envs = InfiniteEnvironments(GLs, GRs)
return recalculate!(envs, above, operator, below; kwargs...)
return recalculate!(envs, below, operator, above; kwargs...)
end

function issamespace(envs::InfiniteEnvironments, above::InfiniteMPS,
function issamespace(envs::InfiniteEnvironments, below::InfiniteMPS,
operator::Union{InfiniteMPO,InfiniteMPOHamiltonian},
below::InfiniteMPS)
L = check_length(above, operator, below)
above::InfiniteMPS)
L = check_length(below, operator, above)
for i in 1:L
space(envs.GLs[i]) ==
(left_virtualspace(below, i) ⊗ left_virtualspace(operator, i)' ←
Expand All @@ -42,39 +42,39 @@ function issamespace(envs::InfiniteEnvironments, above::InfiniteMPS,
return true
end

function recalculate!(envs::InfiniteEnvironments, above::InfiniteMPS,
function recalculate!(envs::InfiniteEnvironments, below::InfiniteMPS,
operator::Union{InfiniteMPO,InfiniteMPOHamiltonian},
below::InfiniteMPS=above; kwargs...)
if !issamespace(envs, above, operator, below)
above::InfiniteMPS=below; kwargs...)
if !issamespace(envs, below, operator, above)
# TODO: in-place initialization?
GLs, GRs = initialize_environments(above, operator, below)
GLs, GRs = initialize_environments(below, operator, above)
copy!(envs.GLs, GLs)
copy!(envs.GRs, GRs)
end

alg = environment_alg(above, operator, below; kwargs...)
alg = environment_alg(below, operator, above; kwargs...)

@sync begin
@spawn compute_leftenvs!(envs, above, operator, below, alg)
@spawn compute_rightenvs!(envs, above, operator, below, alg)
@spawn compute_leftenvs!(envs, below, operator, above, alg)
@spawn compute_rightenvs!(envs, below, operator, above, alg)
end
normalize!(envs, above, operator, below)
normalize!(envs, below, operator, above)

return envs
end

# InfiniteMPO environments
# ------------------------
function initialize_environments(above::InfiniteMPS, operator::InfiniteMPO,
below::InfiniteMPS=above)
L = check_length(above, operator, below)
function initialize_environments(below::InfiniteMPS, operator::InfiniteMPO,
above::InfiniteMPS=below)
L = check_length(below, operator, above)
GLs = PeriodicVector([randomize!(allocate_GL(below, operator, above, i)) for i in 1:L])
GRs = PeriodicVector([randomize!(allocate_GR(below, operator, above, i)) for i in 1:L])
return GLs, GRs
end

function compute_leftenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
operator::InfiniteMPO, below::InfiniteMPS, alg)
function compute_leftenvs!(envs::InfiniteEnvironments, below::InfiniteMPS,
operator::InfiniteMPO, above::InfiniteMPS, alg)
# compute eigenvector
T = TransferMatrix(above.AL, operator, below.AL)
λ, envs.GLs[1] = fixedpoint(flip(T), envs.GLs[1], :LM, alg)
Expand All @@ -86,8 +86,8 @@ function compute_leftenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
return λ, envs
end

function compute_rightenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
operator::InfiniteMPO, below::InfiniteMPS, alg)
function compute_rightenvs!(envs::InfiniteEnvironments, below::InfiniteMPS,
operator::InfiniteMPO, above::InfiniteMPS, alg)
# compute eigenvector
T = TransferMatrix(above.AR, operator, below.AR)
λ, envs.GRs[end] = fixedpoint(T, envs.GRs[end], :LM, alg)
Expand All @@ -99,8 +99,8 @@ function compute_rightenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
return λ, envs
end

function TensorKit.normalize!(envs::InfiniteEnvironments, above::InfiniteMPS,
operator::InfiniteMPO, below::InfiniteMPS)
function TensorKit.normalize!(envs::InfiniteEnvironments, below::InfiniteMPS,
operator::InfiniteMPO, above::InfiniteMPS)
for i in 1:length(operator)
λ = dot(below.C[i], MPO_∂∂C(envs.GLs[i + 1], envs.GRs[i]) * above.C[i])
scale!(envs.GLs[i + 1], inv(λ))
Expand All @@ -110,8 +110,8 @@ end

# InfiniteMPOHamiltonian environments
# -----------------------------------
function initialize_environments(above::InfiniteMPS, operator::InfiniteMPOHamiltonian,
below::InfiniteMPS=above)
function initialize_environments(below::InfiniteMPS, operator::InfiniteMPOHamiltonian,
above::InfiniteMPS=below)
L = check_length(above, operator, below)
GLs = PeriodicVector([allocate_GL(below, operator, above, i) for i in 1:L])
GRs = PeriodicVector([allocate_GR(below, operator, above, i) for i in 1:L])
Expand Down Expand Up @@ -139,9 +139,9 @@ function initialize_environments(above::InfiniteMPS, operator::InfiniteMPOHamilt
return GLs, GRs
end

function compute_leftenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
operator::InfiniteMPOHamiltonian, below::InfiniteMPS, alg)
L = check_length(above, operator, below)
function compute_leftenvs!(envs::InfiniteEnvironments, below::InfiniteMPS,
operator::InfiniteMPOHamiltonian, above::InfiniteMPS, alg)
L = check_length(below, above, operator)
GLs = envs.GLs
vsize = length(first(GLs))

Expand All @@ -156,20 +156,20 @@ function compute_leftenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
# fill_data!(leftutil, one)
# @plansor GL[1][1][-1 -2; -3] = ρ_left[-1; -3] * leftutil[-2]

(L > 1) && left_cyclethrough!(1, GLs, above, operator, below)
(L > 1) && left_cyclethrough!(1, GLs, below, operator, above)

for i in 2:vsize
prev = copy(GLs[1][i])
zerovector!(GLs[1][i])
left_cyclethrough!(i, GLs, above, operator, below)
left_cyclethrough!(i, GLs, below, operator, above)

if isidentitylevel(operator, i) # identity matrices; do the hacky renormalization
T = regularize(TransferMatrix(above.AL, below.AL), ρ_left, ρ_right)
GLs[1][i], convhist = linsolve(flip(T), GLs[1][i], prev, alg, 1, -1)
convhist.converged == 0 &&
@warn "GL$i failed to converge: normres = $(convhist.normres)"

(L > 1) && left_cyclethrough!(i, GLs, above, operator, below)
(L > 1) && left_cyclethrough!(i, GLs, below, operator, above)

# go through the unitcell, again subtracting fixpoints
for site in 1:L
Expand All @@ -186,15 +186,15 @@ function compute_leftenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
convhist.converged == 0 &&
@warn "GL$i failed to converge: normres = $(convhist.normres)"
end
(L > 1) && left_cyclethrough!(i, GLs, above, operator, below)
(L > 1) && left_cyclethrough!(i, GLs, below, operator, above)
end
end

return GLs
end

function left_cyclethrough!(index::Int, GL, above::InfiniteMPS, H::InfiniteMPOHamiltonian,
below::InfiniteMPS=above)
function left_cyclethrough!(index::Int, GL, below::InfiniteMPS, H::InfiniteMPOHamiltonian,
above::InfiniteMPS=below)
# TODO: efficient transfer matrix slicing for large unitcells
leftinds = 1:index
for site in eachindex(GL)
Expand All @@ -205,8 +205,8 @@ function left_cyclethrough!(index::Int, GL, above::InfiniteMPS, H::InfiniteMPOHa
return GL
end

function compute_rightenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
operator::InfiniteMPOHamiltonian, below::InfiniteMPS, alg)
function compute_rightenvs!(envs::InfiniteEnvironments, below::InfiniteMPS,
operator::InfiniteMPOHamiltonian, above::InfiniteMPS, alg)
L = check_length(above, operator, below)
GRs = envs.GRs
vsize = length(last(GRs))
Expand All @@ -222,12 +222,12 @@ function compute_rightenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
# fill_data!(rightutil, one)
# @plansor GR[end][end][-1 -2; -3] = r_RR(state)[-1; -3] * rightutil[-2]

(L > 1) && right_cyclethrough!(vsize, GRs, above, operator, below) # populate other sites
(L > 1) && right_cyclethrough!(vsize, GRs, below, operator, above) # populate other sites

for i in (vsize - 1):-1:1
prev = copy(GRs[end][i])
zerovector!(GRs[end][i])
right_cyclethrough!(i, GRs, above, operator, below)
right_cyclethrough!(i, GRs, below, operator, above)

if isidentitylevel(operator, i) # identity matrices; do the hacky renormalization
# subtract fixpoints
Expand All @@ -236,7 +236,7 @@ function compute_rightenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
convhist.converged == 0 &&
@warn "GR$i failed to converge: normres = $(convhist.normres)"

L > 1 && right_cyclethrough!(i, GRs, above, operator, below)
L > 1 && right_cyclethrough!(i, GRs, below, operator, above)

# go through the unitcell, again subtracting fixpoints
for site in 1:L
Expand All @@ -253,16 +253,16 @@ function compute_rightenvs!(envs::InfiniteEnvironments, above::InfiniteMPS,
@warn "GR$i failed to converge: normres = $(convhist.normres)"
end

(L > 1) && right_cyclethrough!(i, GRs, above, operator, below)
(L > 1) && right_cyclethrough!(i, GRs, below, operator, above)
end
end

return GRs
end

function right_cyclethrough!(index::Int, GR, above::InfiniteMPS,
function right_cyclethrough!(index::Int, GR, below::InfiniteMPS,
operator::InfiniteMPOHamiltonian,
below::InfiniteMPS)
above::InfiniteMPS=below)
# TODO: efficient transfer matrix slicing for large unitcells
for site in reverse(eachindex(GR))
rightinds = index:length(GR[site])
Expand All @@ -274,21 +274,21 @@ function right_cyclethrough!(index::Int, GR, above::InfiniteMPS,
end

# no normalization necessary -- for consistant interface
function TensorKit.normalize!(envs::InfiniteEnvironments, above::InfiniteMPS,
operator::InfiniteMPOHamiltonian, below::InfiniteMPS)
function TensorKit.normalize!(envs::InfiniteEnvironments, below::InfiniteMPS,
operator::InfiniteMPOHamiltonian, above::InfiniteMPS)
return envs
end

# Transfer operations
# -------------------

function transfer_leftenv!(envs::InfiniteEnvironments, above, operator, below, site::Int)
function transfer_leftenv!(envs::InfiniteEnvironments, below, operator, above, site::Int)
T = TransferMatrix(above.AL[site - 1], operator[site - 1], below.AL[site - 1])
envs.GLs[site] = envs.GLs[site - 1] * T
return envs
end

function transfer_rightenv!(envs::InfiniteEnvironments, above, operator, below, site::Int)
function transfer_rightenv!(envs::InfiniteEnvironments, below, operator, above, site::Int)
T = TransferMatrix(above.AR[site + 1], operator[site + 1], below.AR[site + 1])
envs.GRs[site] = T * envs.GRs[site + 1]
return envs
Expand Down
34 changes: 17 additions & 17 deletions src/environments/multiline_envs.jl
Original file line number Diff line number Diff line change
@@ -1,38 +1,38 @@
const MultilineEnvironments{E<:AbstractMPSEnvironments} = Multiline{E}

function environments(above::MultilineMPS, operator::MultilineMPO,
below::MultilineMPS=above; kwargs...)
function environments(below::MultilineMPS, operator::MultilineMPO,
above::MultilineMPS=below; kwargs...)
(rows = size(above, 1)) == size(operator, 1) == size(below, 1) ||
throw(ArgumentError("Incompatible sizes"))
envs = map(1:rows) do row
return environments(above[row], operator[row], below[row + 1]; kwargs...)
return environments(below[row + 1], operator[row], above[row]; kwargs...)
end
return Multiline(PeriodicVector(envs))
end

function recalculate!(envs::MultilineEnvironments, above::MultilineMPS,
operator::MultilineMPO, below::MultilineMPS=above; kwargs...)
function recalculate!(envs::MultilineEnvironments, below::MultilineMPS,
operator::MultilineMPO, above::MultilineMPS=below; kwargs...)
(rows = size(above, 1)) == size(operator, 1) == size(below, 1) ||
throw(ArgumentError("Incompatible sizes"))
@threads for row in 1:rows
recalculate!(envs[row], above[row], operator[row], below[row + 1]; kwargs...)
recalculate!(envs[row], below[row + 1], operator[row], above[row]; kwargs...)
end
return envs
end
function recalculate!(envs::MultilineEnvironments, below, (operator, above)::Tuple;
kwargs...)
return recalculate!(envs, above, operator, below; kwargs...)
return recalculate!(envs, below, operator, above; kwargs...)
end

function TensorKit.normalize!(envs::MultilineEnvironments, above, operator, below)
for row in 1:size(above, 1)
normalize!(envs[row], above[row], operator[row], below[row + 1])
function TensorKit.normalize!(envs::MultilineEnvironments, below, operator, above)
for row in 1:size(below, 1)
normalize!(envs[row], below[row + 1], operator[row], above[row])
end
return envs
end
function TensorKit.normalize!(envs::MultilineEnvironments, below, (operator, above))
for row in 1:size(above, 1)
normalize!(envs[row], above[row], operator[row], below[row + 1])
normalize!(envs[row], below[row + 1], operator[row], above[row])
end
return envs
end
Expand All @@ -44,22 +44,22 @@ function rightenv(envs::MultilineEnvironments, col::Int, state)
return rightenv.(parent(envs), col, parent(state))
end

function transfer_leftenv!(envs::MultilineEnvironments, above, operator, below, site::Int)
function transfer_leftenv!(envs::MultilineEnvironments, below, operator, above, site::Int)
for row in 1:size(above, 1)
transfer_leftenv!(envs[row], above[row], operator[row], below[row + 1], site)
transfer_leftenv!(envs[row], below[row + 1], operator[row], above[row], site)
end
return envs
end
function transfer_leftenv!(envs::MultilineEnvironments, below, (O, above)::Tuple, site::Int)
return transfer_leftenv!(envs, above, O, below, site)
return transfer_leftenv!(envs, below, O, above, site)
end
function transfer_rightenv!(envs::MultilineEnvironments, above, operator, below, site::Int)
function transfer_rightenv!(envs::MultilineEnvironments, below, operator, above, site::Int)
for row in 1:size(above, 1)
transfer_rightenv!(envs[row], above[row], operator[row], below[row + 1], site)
transfer_rightenv!(envs[row], below[row + 1], operator[row], above[row], site)
end
return envs
end
function transfer_rightenv!(envs::MultilineEnvironments, below, (O, above)::Tuple,
site::Int)
return transfer_rightenv!(envs, above, O, below, site)
return transfer_rightenv!(envs, below, O, above, site)
end
Loading