Update for MPSKit v0.12 (#42)

* Update for BlockTensors

* Bugfix mpoham constructor

adds missing end of each term

* update compat

* `AbstractTensorMap` typesignature changes

* always use `AbstractTensorMap` as tensortype of `MPOHamiltonian`s, correct braidingtensor signature

* use `SparseBlockTensorMap`s for `MPOHamiltonian`s

* allow for different types of tensors at different sites

* add new `AbstractFiniteLattice` and `AbstractInfiniteLattice` supertypes for differentiating between finite and infinite models

* implement `Base.isfinite` for lattice(types)

* add `Base.isfinite` for LatticePoints

* make sure finite lattices result in a `FiniteMPOHamiltonian`

* change the handling of checking whether a lattice is finite, remove abstract types introduced earlier

* account for changes to `right_virtualspace`

* update TensorKit compat

* Bump version v0.4.0

* Small fixes

* Update spinoperators for TensorKit version

* update bosonoperators for new TensorKit version

* update fermionoperators for new TensorKit version

* update hubbardoperators for TensorKit version

* update tJ operators for TensorKit version

* Revert "allow for different types of tensors at different sites"

This reverts commit a09994c.

* rework mpoham to use existing functionality

* disable type inference testing

* Cleanup

* update CI

---------

Co-authored-by: victor <[email protected]>

Author: lkdvos

.github/workflows/CI.yml

@@ -21,20 +21,14 @@ jobs:
       fail-fast: false
       matrix:
         version:
-          - '1.8' # LTS version
+          - 'lts' # LTS version
           - '1' # automatically expands to the latest stable 1.x release of Julia
         os:
           - ubuntu-latest
           - macOS-latest
           - windows-latest
         arch:
           - x64
-          # - x86 # exclude x86 tests - not enough memory and already failing on TensorKt
-        # exclude:
-        #   - os: macOS-latest
-        #     arch: x86
-        #   - os: windows-latest
-        #     arch: x86
     steps:
       - uses: actions/checkout@v4
       - uses: julia-actions/setup-julia@v2
@@ -52,27 +46,3 @@ jobs:
           CODECOV_TOKEN: ${{ secrets.CODECOV_TOKEN }}
         with:
           file: lcov.info
-  test-nightly:
-    needs: test
-    name: Julia nightly - ${{ matrix.os }} - ${{ matrix.arch }}
-    runs-on: ${{ matrix.os }}
-    strategy:
-      fail-fast: false
-      matrix:
-        version:
-          - 'nightly'
-        os:
-          - ubuntu-latest
-        arch:
-          - x64
-    steps:
-      - uses: actions/checkout@v4
-      - uses: julia-actions/setup-julia@v2
-        with:
-          version: ${{ matrix.version }}
-          arch: ${{ matrix.arch }}
-      - uses: julia-actions/cache@v2
-      - uses: julia-actions/julia-buildpkg@latest
-      - uses: julia-actions/julia-runtest@latest
-        env:
-          JULIA_NUM_THREADS: 4

Project.toml

@@ -1,7 +1,7 @@
 name = "MPSKitModels"
 uuid = "ca635005-6f8c-4cd1-b51d-8491250ef2ab"
 authors = ["Maarten Van Damme", "Lukas Devos", "Gertian Roose", "Klaas Gunst"]
-version = "0.3.6"
+version = "0.4.0"
 
 [deps]
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
@@ -13,13 +13,13 @@ TensorOperations = "6aa20fa7-93e2-5fca-9bc0-fbd0db3c71a2"
 TupleTools = "9d95972d-f1c8-5527-a6e0-b4b365fa01f6"
 
 [compat]
-MPSKit = "0.11"
+MPSKit = "0.12"
 MacroTools = "0.5"
 PrecompileTools = "1"
-TensorKit = "0.11, 0.12"
-TensorOperations = "4"
+TensorKit = "0.13,0.14"
+TensorOperations = "5"
 TupleTools = "1"
-julia = "1.8"
+julia = "1.10"
 
 [extras]
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

src/MPSKitModels.jl

@@ -2,7 +2,7 @@ module MPSKitModels
 
 using TensorKit, MPSKit
 using MacroTools: @capture, postwalk
-using MPSKit: @plansor, _lastspace, _firstspace
+using MPSKit: @plansor
 using TensorOperations
 using TupleTools
 

src/lattices/chains.jl

@@ -25,6 +25,7 @@ struct FiniteChain <: AbstractLattice{1}
 end
 Base.axes(chain::FiniteChain) = (1:(chain.L),)
 Base.isfinite(::Type{FiniteChain}) = true
+Base.isfinite(::FiniteChain) = true
 
 const Chain = Union{InfiniteChain,FiniteChain}
 

src/lattices/latticepoints.jl

@@ -61,5 +61,7 @@ end
 Base.:-(i::Int, j::LatticePoint{1}) = LatticePoint(i .- j, j.lattice)
 
 Base.isless(i::L, j::L) where {L<:LatticePoint} = linearize_index(i) < linearize_index(j)
-
+function Base.isfinite(::Union{LatticePoint{N,G},Type{<:LatticePoint{N,G}}}) where {N,G}
+    return isfinite(G)
+end
 latticetype(::Union{LatticePoint{N,G},Type{<:LatticePoint{N,G}}}) where {N,G} = G

src/models/hamiltonians.jl

@@ -40,8 +40,8 @@ function transverse_field_ising(T::Type{<:Number}=ComplexF64,
                                 S::Union{Type{Trivial},Type{Z2Irrep}}=Trivial,
                                 lattice::AbstractLattice=InfiniteChain(1);
                                 J=1.0, g=1.0)
-    ZZ = rmul!(σᶻᶻ(T, S), -J)
-    X = rmul!(σˣ(T, S), g * -J)
+    ZZ = scale!(σᶻᶻ(T, S), -J)
+    X = scale!(σˣ(T, S), g * -J)
     return @mpoham begin
         sum(nearest_neighbours(lattice)) do (i, j)
             return ZZ{i,j}
@@ -53,8 +53,13 @@ end
 function transverse_field_ising(T::Type{<:Number}, ::Type{fℤ₂},
                                 lattice::AbstractLattice=InfiniteChain(1);
                                 J=1.0, g=1.0)
-    twosite = axpby!(-J, c_plusmin(T) + c_minplus(T), J, c_plusplus(T) + c_minmin(T))
-    onesite = axpby!(2g * J, c_number(T), -g * J, id(Matrix{T}, space(twosite, 1)))
+    hop = add!(c_plusmin(T), c_minplus(T))
+    sc = add!(c_plusplus(T), c_minmin(T))
+    twosite = add!(hop, sc, J, -J)
+
+    N = c_number(T)
+    E = id(storagetype(N), space(N, 1))
+    onesite = add!(N, E, -g * J, 2g * J)
 
     return @mpoham begin
         sum(nearest_neighbours(lattice)) do (i, j)

src/models/transfermatrices.jl

@@ -30,19 +30,19 @@ function classical_ising(elt::Type{<:Number}=ComplexF64, ::Type{Trivial}=Trivial
 
     @tensor o[-1 -2; -3 -4] := O[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * nt[-3; 3] * nt[-4; 4]
 
-    return DenseMPO(TensorMap(o, ℂ^2 * ℂ^2, ℂ^2 * ℂ^2))
+    return InfiniteMPO([TensorMap(o, ℂ^2 * ℂ^2, ℂ^2 * ℂ^2)])
 end
 
 function classical_ising(elt::Type{<:Number}, ::Type{Z2Irrep}; beta=log(1 + sqrt(2)) / 2)
     x = cosh(beta)
     y = sinh(beta)
 
     sec = ℤ₂Space(0 => 1, 1 => 1)
-    mpo = TensorMap(zeros, elt, sec * sec, sec * sec)
-    blocks(mpo)[Irrep[ℤ₂](0)] = [2x^2 2x*y; 2x*y 2y^2]
-    blocks(mpo)[Irrep[ℤ₂](1)] = [2x*y 2x*y; 2x*y 2x*y]
+    mpo = zeros(elt, sec * sec, sec * sec)
+    block(mpo, Irrep[ℤ₂](0)) .= [2x^2 2x*y; 2x*y 2y^2]
+    block(mpo, Irrep[ℤ₂](1)) .= [2x*y 2x*y; 2x*y 2x*y]
 
-    return DenseMPO(mpo)
+    return InfiniteMPO([mpo])
 end
 
 #===========================================================================================
@@ -63,23 +63,23 @@ function sixvertex(elt::Type{<:Number}=ComplexF64, ::Type{Trivial}=Trivial; a=1.
             0 c b 0
             0 b c 0
             0 0 0 a]
-    return DenseMPO(permute(TensorMap(d, ℂ^2 ⊗ ℂ^2, ℂ^2 ⊗ ℂ^2), ((1, 2), (4, 3))))
+    return InfiniteMPO([permute(TensorMap(d, ℂ^2 ⊗ ℂ^2, ℂ^2 ⊗ ℂ^2), ((1, 2), (4, 3)))])
 end
 function sixvertex(elt::Type{<:Number}, ::Type{U1Irrep}; a=1.0, b=1.0, c=1.0)
     pspace = U1Space(-1 // 2 => 1, 1 // 2 => 1)
-    mpo = TensorMap(zeros, elt, pspace ⊗ pspace, pspace ⊗ pspace)
-    blocks(mpo)[Irrep[U₁](0)] = [b c; c b]
-    blocks(mpo)[Irrep[U₁](1)] = reshape([a], (1, 1))
-    blocks(mpo)[Irrep[U₁](-1)] = reshape([a], (1, 1))
-    return DenseMPO(permute(mpo, ((1, 2), (4, 3))))
+    mpo = zeros(elt, pspace ⊗ pspace, pspace ⊗ pspace)
+    block(mpo, Irrep[U₁](0)) .= [b c; c b]
+    block(mpo, Irrep[U₁](1)) .= reshape([a], (1, 1))
+    block(mpo, Irrep[U₁](-1)) .= reshape([a], (1, 1))
+    return InfiniteMPO([permute(mpo, ((1, 2), (4, 3)))])
 end
 function sixvertex(elt::Type{<:Number}, ::Type{CU1Irrep}; a=1.0, b=1.0, c=1.0)
     pspace = CU1Space(1 // 2 => 1)
-    mpo = TensorMap(zeros, elt, pspace ⊗ pspace, pspace ⊗ pspace)
-    blocks(mpo)[Irrep[CU₁](0, 0)] = reshape([b + c], (1, 1))
-    blocks(mpo)[Irrep[CU₁](0, 1)] = reshape([-b + c], (1, 1))
-    blocks(mpo)[Irrep[CU₁](1, 2)] = reshape([a], (1, 1))
-    return DenseMPO(permute(mpo, ((1, 2), (4, 3))))
+    mpo = zeros(elt, pspace ⊗ pspace, pspace ⊗ pspace)
+    block(mpo, Irrep[CU₁](0, 0)) .= reshape([b + c], (1, 1))
+    block(mpo, Irrep[CU₁](0, 1)) .= reshape([-b + c], (1, 1))
+    block(mpo, Irrep[CU₁](1, 2)) .= reshape([a], (1, 1))
+    return InfiniteMPO([permute(mpo, ((1, 2), (4, 3)))])
 end
 
 #===========================================================================================
@@ -93,9 +93,9 @@ MPO for the partition function of the two-dimensional hard hexagon model.
 """
 function hard_hexagon(elt::Type{<:Number}=ComplexF64)
     P = Vect[FibonacciAnyon](:τ => 1)
-    O = TensorMap(ones, elt, P ⊗ P ← P ⊗ P)
-    blocks(O)[FibonacciAnyon(:I)] *= 0
-    return DenseMPO(O)
+    O = ones(elt, P ⊗ P ← P ⊗ P)
+    block(O, FibonacciAnyon(:I)) .*= 0
+    return InfiniteMPO([O])
 end
 
 #===========================================================================================
@@ -116,5 +116,5 @@ function qstate_clock(elt::Type{<:Number}=ComplexF64, ::Type{Trivial}=Trivial;
                             (comega(i - j) + comega(j - k) + comega(k - l) + comega(l - i)))
     end
 
-    return DenseMPO(TensorMap(O, ℂ^q * ℂ^q, ℂ^q * ℂ^q))
+    return InfiniteMPO([TensorMap(O, ℂ^q * ℂ^q, ℂ^q * ℂ^q)])
 end

src/operators/bosonoperators.jl

@@ -10,7 +10,7 @@ a_plus(elt::Type{<:Number}; kwargs...) = a_plus(elt, Trivial; kwargs...)
 a_plus(symm::Type{<:Sector}; kwargs...) = a_plus(ComplexF64, symm; kwargs...)
 
 function a_plus(elt::Type{<:Number}, ::Type{Trivial}; cutoff::Integer=5)
-    a⁺ = TensorMap(zeros, elt, ComplexSpace(cutoff + 1), ComplexSpace(cutoff + 1))
+    a⁺ = zeros(elt, ComplexSpace(cutoff + 1), ComplexSpace(cutoff + 1))
     for n in 1:cutoff
         a⁺[n + 1, n] = sqrt(n)
     end
@@ -21,7 +21,7 @@ function a_plus(elt::Type{<:Number}, ::Type{U1Irrep}; cutoff::Integer=5, side=:L
     pspace = U1Space(n => 1 for n in 0:cutoff)
     if side === :L
         vspace = U1Space(1 => 1)
-        a⁺ = TensorMap(zeros, elt, pspace ← pspace ⊗ vspace)
+        a⁺ = zeros(elt, pspace ← pspace ⊗ vspace)
         for (f1, f2) in fusiontrees(a⁺)
             c₁, c₂ = f1.uncoupled[1], f2.uncoupled[1]
             if c₁.charge == c₂.charge + 1
@@ -30,7 +30,7 @@ function a_plus(elt::Type{<:Number}, ::Type{U1Irrep}; cutoff::Integer=5, side=:L
         end
     elseif side === :R
         vspace = U1Space(-1 => 1)
-        a⁺ = TensorMap(zeros, elt, vspace ⊗ pspace ← pspace)
+        a⁺ = zeros(elt, vspace ⊗ pspace ← pspace)
         for (f1, f2) in fusiontrees(a⁺)
             c₁, c₂ = f1.uncoupled[2], f2.uncoupled[1]
             if c₁.charge == c₂.charge + 1
@@ -57,7 +57,7 @@ a_min(elt::Type{<:Number}; kwargs...) = a_min(elt, Trivial; kwargs...)
 a_min(symm::Type{<:Sector}; kwargs...) = a_min(ComplexF64, symm; kwargs...)
 
 function a_min(elt::Type{<:Number}, ::Type{Trivial}; cutoff::Integer=5)
-    a⁻ = TensorMap(zeros, elt, ComplexSpace(cutoff + 1), ComplexSpace(cutoff + 1))
+    a⁻ = zeros(elt, ComplexSpace(cutoff + 1), ComplexSpace(cutoff + 1))
     for n in 1:cutoff
         a⁻[n, n + 1] = sqrt(n)
     end
@@ -68,7 +68,7 @@ function a_min(elt::Type{<:Number}, ::Type{U1Irrep}; cutoff::Integer=5, side=:L)
     pspace = U1Space(n => 1 for n in 0:cutoff)
     if side === :L
         vspace = U1Space(-1 => 1)
-        a⁻ = TensorMap(zeros, elt, pspace ← pspace ⊗ vspace)
+        a⁻ = zeros(elt, pspace ← pspace ⊗ vspace)
         for (f1, f2) in fusiontrees(a⁻)
             c₁, c₂ = f1.uncoupled[1], f2.uncoupled[1]
             if c₁.charge + 1 == c₂.charge
@@ -77,7 +77,7 @@ function a_min(elt::Type{<:Number}, ::Type{U1Irrep}; cutoff::Integer=5, side=:L)
         end
     elseif side === :R
         vspace = U1Space(1 => 1)
-        a⁻ = TensorMap(zeros, elt, vspace ⊗ pspace ← pspace)
+        a⁻ = zeros(elt, vspace ⊗ pspace ← pspace)
         for (f1, f2) in fusiontrees(a⁻)
             c₁, c₂ = f1.uncoupled[2], f2.uncoupled[1]
             if c₁.charge + 1 == c₂.charge
@@ -123,7 +123,7 @@ a_number(elt::Type{<:Number}; kwargs...) = a_number(elt, Trivial; kwargs...)
 a_number(symm::Type{<:Sector}; kwargs...) = a_number(ComplexF64, symm; kwargs...)
 
 function a_number(elt::Type{<:Number}, ::Type{Trivial}; cutoff::Integer=5)
-    N = TensorMap(zeros, elt, ComplexSpace(cutoff + 1), ComplexSpace(cutoff + 1))
+    N = zeros(elt, ComplexSpace(cutoff + 1), ComplexSpace(cutoff + 1))
     for n in 0:cutoff
         N[n + 1, n + 1] = n
     end
@@ -132,7 +132,7 @@ end
 
 function a_number(elt::Type{<:Number}, ::Type{U1Irrep}; cutoff::Integer=5)
     pspace = U1Space(n => 1 for n in 0:cutoff)
-    N = TensorMap(zeros, elt, pspace, pspace)
+    N = zeros(elt, pspace, pspace)
     for (c, b) in blocks(N)
         b .= c.charge
     end

src/operators/fermionoperators.jl

@@ -12,8 +12,8 @@ function c_plus(elt::Type{<:Number}=ComplexF64; side=:L)
     vspace = Vect[fℤ₂](1 => 1)
     if side === :L
         pspace = Vect[fℤ₂](0 => 1, 1 => 1)
-        c⁺ = TensorMap(zeros, elt, pspace ← pspace ⊗ vspace)
-        blocks(c⁺)[fℤ₂(1)] .= one(elt)
+        c⁺ = zeros(elt, pspace ← pspace ⊗ vspace)
+        block(c⁺, fℤ₂(1)) .= one(elt)
     elseif side === :R
         C = c_plus(elt; side=:L)
         F = isomorphism(storagetype(C), vspace, flip(vspace))
@@ -62,7 +62,7 @@ Fermionic number operator.
 """
 function c_number(elt::Type{<:Number}=ComplexF64)
     pspace = Vect[fℤ₂](0 => 1, 1 => 1)
-    n = TensorMap(zeros, elt, pspace ← pspace)
-    blocks(n)[fℤ₂(1)] .= one(elt)
+    n = zeros(elt, pspace ← pspace)
+    block(n, fℤ₂(1)) .= one(elt)
     return n
 end

src/operators/hubbardoperators.jl

@@ -52,13 +52,13 @@ end
 function single_site_operator(T, particle_symmetry::Type{<:Sector},
                               spin_symmetry::Type{<:Sector})
     V = hubbard_space(particle_symmetry, spin_symmetry)
-    return TensorMap(zeros, T, V ← V)
+    return zeros(T, V ← V)
 end
 
 function two_site_operator(T, particle_symmetry::Type{<:Sector},
                            spin_symmetry::Type{<:Sector})
     V = hubbard_space(particle_symmetry, spin_symmetry)
-    return TensorMap(zeros, T, V ⊗ V ← V ⊗ V)
+    return zeros(T, V ⊗ V ← V ⊗ V)
 end
 
 """