windowhamiltonian poc

src/MPSKit.jl

@@ -19,7 +19,7 @@ export QP, LeftGaugedQP, RightGaugedQP
 export AbstractMPO
 export MPO, FiniteMPO, InfiniteMPO
 export JordanMPOTensor, JordanMPOTensorMap
-export MPOHamiltonian, FiniteMPOHamiltonian, InfiniteMPOHamiltonian
+export MPOHamiltonian, FiniteMPOHamiltonian, InfiniteMPOHamiltonian, WindowMPOHamiltonian
 export MultilineMPO
 export UntimedOperator, TimedOperator, MultipliedOperator, LazySum
 
@@ -122,6 +122,7 @@ include("operators/abstractmpo.jl")
 include("operators/mpo.jl")
 include("operators/jordanmpotensor.jl")
 include("operators/mpohamiltonian.jl") # the mpohamiltonian objects
+include("operators/windowhamiltonian.jl")
 include("operators/ortho.jl")
 include("operators/multilinempo.jl")
 include("operators/projection.jl")

src/algorithms/expval.jl

@@ -128,6 +128,13 @@ function expectation_value(
     return dot(ψ, H, ψ, envs) / dot(ψ, ψ)
 end
 
+function expectation_value(
+        ψ::WindowMPS, H::WindowMPOHamiltonian,
+        envs::AbstractMPSEnvironments = environments(ψ, H)
+    )
+    return dot(ψ, H, ψ, envs) / dot(ψ, ψ)
+end
+
 function expectation_value(
         ψ::InfiniteMPS, H::InfiniteMPOHamiltonian,
         envs::AbstractMPSEnvironments = environments(ψ, H)

src/algorithms/propagator/corvector.jl

@@ -59,7 +59,7 @@ See also [`Jeckelmann`](@ref) for the original approach.
 struct NaiveInvert <: DDMRG_Flavour end
 
 function propagator(
-        A::AbstractFiniteMPS, z::Number, H::FiniteMPOHamiltonian,
+        A::AbstractFiniteMPS, z::Number, H,
         alg::DynamicalDMRG{NaiveInvert}; init = copy(A)
     )
     h_envs = environments(init, H) # environments for h
@@ -119,7 +119,7 @@ See also [`NaiveInvert`](@ref) for a less costly but less accurate alternative.
 struct Jeckelmann <: DDMRG_Flavour end
 
 function propagator(
-        A::AbstractFiniteMPS, z, H::FiniteMPOHamiltonian,
+        A::AbstractFiniteMPS, z, H,
         alg::DynamicalDMRG{Jeckelmann}; init = copy(A)
     )
     ω = real(z)
@@ -176,7 +176,7 @@ function propagator(
 end
 
 function squaredenvs(
-        state::AbstractFiniteMPS, H::FiniteMPOHamiltonian, envs = environments(state, H)
+        state::AbstractFiniteMPS, H, envs = environments(state, H)
     )
     H² = conj(H) * H
     L = length(state)

src/environments/finite_envs.jl

@@ -53,14 +53,14 @@ function environments(
     return environments(below, O, above, leftstart, rightstart)
 end
 function environments(
-        below::WindowMPS, O::Union{InfiniteMPOHamiltonian, InfiniteMPO}, above = nothing;
-        lenvs = environments(below.left_gs, O),
-        renvs = environments(below.right_gs, O)
+        below::WindowMPS, O::WindowMPOHamiltonian, above = nothing;
+        lenvs = environments(below.left_gs, O.left_ham),
+        renvs = environments(below.right_gs, O.right_ham)
     )
     leftstart = copy(lenvs.GLs[1])
     rightstart = copy(renvs.GRs[end])
 
-    return environments(below, O, above, leftstart, rightstart)
+    return environments(below, O.finite_ham, above, leftstart, rightstart)
 end
 
 function environments(below::S, above::S) where {S <: Union{FiniteMPS, WindowMPS}}

src/operators/windowhamiltonian.jl

@@ -0,0 +1,40 @@
+"""
+A WindowMPS is a finite MPS embedded between an infinite mps to the left, and an inifite mps to the right. 
+The associated hamiltonian has also an infinite part to the left, a finite hamiltonian in the middle, and an infinite part to the right.
+
+Acts simalar as just a finite hamiltonian, but we 'remember' the boundary hamiltonians.
+"""
+
+# todo - what is the required interface for abstractmpo?
+# support densempo windows?
+struct WindowMPOHamiltonian{O} <: AbstractMPO{O}
+    left_ham::InfiniteMPOHamiltonian{O}
+    finite_ham::FiniteMPOHamiltonian{O}
+    right_ham::InfiniteMPOHamiltonian{O}
+end
+
+#utility constructor
+function WindowMPOHamiltonian(ham::InfiniteMPOHamiltonian, interval::UnitRange)
+
+    # to make sure the interval corresponds with finite_ham, it is important that the unitcell of the left/right hamiltonians is circshifted correctly
+    left_edge = (interval.start - 1) % length(ham)
+    left_ham = InfiniteMPOHamiltonian([ham[i] for i in (left_edge - length(ham) + 1):left_edge])
+    right_edge = (interval.stop + 1) % length(ham)
+    right_ham = InfiniteMPOHamiltonian([ham[i] for i in right_edge:(right_edge + length(ham) - 1)])
+
+    finite_ham = FiniteMPOHamiltonian([ham[i] for i in interval])
+    return WindowMPOHamiltonian(left_ham, finite_ham, right_ham)
+end
+
+
+Base.copy(h::WindowMPOHamiltonian) = WindowMPOHamiltonian(copy(h.left_ham), copy(h.finite_ham), copy(h.right_ham))
+
+# some basic linalg
+for fun in (:(Base.:+), :(Base.:-), :(Base.:*))
+    @eval $fun(a::WindowMPOHamiltonian, b::WindowMPOHamiltonian) = WindowMPOHamiltonian($fun(a.left_ham, b.left_ham), $fun(a.finite_ham, b.finite_ham), $fun(a.right_ham, b.right_ham))
+end
+
+TensorKit.dot(
+    bra::WindowMPS, H::WindowMPOHamiltonian, ket::WindowMPS = bra,
+    envs = environments(bra, H, ket)
+) = dot(bra.window, H.finite_ham, ket.window, envs)

test/algorithms/dynamical_dmrg.jl

@@ -12,7 +12,7 @@ using TensorKit: ℙ
 
 verbosity_conv = 1
 
-@testset "Dynamical DMRG" verbose = true begin
+@testset "Dynamical DMRG FiniteMPS" verbose = true begin
     L = 10
     H = force_planar(-transverse_field_ising(; L, g = -4))
     gs, = find_groundstate(FiniteMPS(L, ℙ^2, ℙ^10), H; verbosity = verbosity_conv)
@@ -32,3 +32,32 @@ verbosity_conv = 1
         @test data ≈ predicted atol = 1.0e-8
     end
 end
+
+
+@testset "Dynamical DMRG WindowMPS" verbose = true begin
+    N = 20
+
+    H = transverse_field_ising(g = -4)
+    Ω = InfiniteMPS(ComplexSpace(2), ComplexSpace(20))
+
+    (Ω, _) = find_groundstate(Ω, H, VUMPS(verbosity = verbosity_conv))
+    XΩ = WindowMPS(Ω, N)
+    H_w = WindowMPOHamiltonian(H, 1:N)
+
+
+    gs_en = expectation_value(XΩ, H_w)
+    vals = range(gs_en - 1.0, gs_en + 1.0, length = 5)
+    eta = 0.3im
+    predicted = [1 / (v + eta - gs_en) for v in vals]
+
+
+    @testset "Flavour $f" for f in (Jeckelmann(), NaiveInvert())
+        alg = DynamicalDMRG(; flavour = f, verbosity = 0, tol = 1.0e-8)
+        data = map(vals) do v
+            result, = propagator(XΩ, v + eta, H_w, alg)
+            return result
+        end
+        @test data ≈ predicted atol = 1.0e-8
+    end
+end
+

test/states/windowmps.jl

@@ -63,14 +63,15 @@ using TensorKit: ℙ
 
     e1 = expectation_value(window, (2, 3) => O)
 
-    window, envs, _ = find_groundstate(window, ham, DMRG(; verbosity = 0))
+    w_ham = WindowMPOHamiltonian(ham, 1:length(window))
+    window, envs, _ = find_groundstate(window, w_ham, DMRG(; verbosity = 0))
 
     e2 = expectation_value(window, (2, 3) => O)
 
     @test real(e2) ≤ real(e1)
 
-    window, envs = timestep(window, ham, 0.1, 0.0, TDVP2(; trscheme = truncrank(20)), envs)
-    window, envs = timestep(window, ham, 0.1, 0.0, TDVP(), envs)
+    window, envs = timestep(window, w_ham, 0.1, 0.0, TDVP2(; trscheme = truncrank(20)), envs)
+    window, envs = timestep(window, w_ham, 0.1, 0.0, TDVP(), envs)
 
     e3 = expectation_value(window, (2, 3) => O)