Merge branch 'main' into ITensors_fix_eigen_perturbation

Author: mtfishman

ITensorGPU/test/test_dmrg.jl

@@ -93,7 +93,8 @@ end
     dmrg_params = (; nsweeps=4, maxdim=10, cutoff=1e-8, noise=1e-8, outputlevel=0)
     observer = DMRGObserver(["Z"], s; energy_tol=1e-4, minsweeps=10)
     E, ψ = dmrg(H, ψ0; observer=observer, dmrg_params...)
-    @test expect(ψ, "Z") ≈ observer.measurements["Z"][end]
+    @test expect(ψ, "Z") ≈ observer.measurements["Z"][end] rtol =
+      10 * sqrt(eps(real(ITensors.scalartype(ψ0))))
     @test correlation_matrix(ψ, "Z", "Z") ≈ correlation_matrix(cpu(ψ), "Z", "Z")
   end
 

NDTensors/Project.toml

@@ -12,10 +12,12 @@ Compat = "34da2185-b29b-5c13-b0c7-acf172513d20"
 Dictionaries = "85a47980-9c8c-11e8-2b9f-f7ca1fa99fb4"
 EllipsisNotation = "da5c29d0-fa7d-589e-88eb-ea29b0a81949"
 FLoops = "cc61a311-1640-44b5-9fba-1b764f453329"
+FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b"
 Folds = "41a02a25-b8f0-4f67-bc48-60067656b558"
 Functors = "d9f16b24-f501-4c13-a1f2-28368ffc5196"
 GPUArraysCore = "46192b85-c4d5-4398-a991-12ede77f4527"
 HDF5 = "f67ccb44-e63f-5c2f-98bd-6dc0ccc4ba2f"
+HalfIntegers = "f0d1745a-41c9-11e9-1dd9-e5d34d218721"
 InlineStrings = "842dd82b-1e85-43dc-bf29-5d0ee9dffc48"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 MappedArrays = "dbb5928d-eab1-5f90-85c2-b9b0edb7c900"
@@ -52,10 +54,12 @@ Compat = "4.9"
 Dictionaries = "0.3.5"
 EllipsisNotation = "1.8"
 FLoops = "0.2.1"
+FillArrays = "1"
 Folds = "0.2.8"
 Functors = "0.2, 0.3, 0.4"
 GPUArraysCore = "0.1"
 HDF5 = "0.14, 0.15, 0.16, 0.17"
+HalfIntegers = "1"
 InlineStrings = "1"
 LinearAlgebra = "1.6"
 MappedArrays = "0.4"

NDTensors/src/NDTensors.jl

@@ -23,6 +23,7 @@ for lib in [
   :AlgorithmSelection,
   :AllocateData,
   :BaseExtensions,
+  :UnspecifiedTypes,
   :SetParameters,
   :BroadcastMapConversion,
   :Unwrap,
@@ -35,8 +36,10 @@ for lib in [
   :BlockSparseArrays,
   :NamedDimsArrays,
   :SmallVectors,
+  :Sectors,
   :SortedSets,
   :TagSets,
+  :UnallocatedArrays,
 ]
   include("lib/$(lib)/src/$(lib).jl")
   @eval using .$lib: $lib

NDTensors/src/lib/Sectors/src/Sectors.jl

@@ -0,0 +1,10 @@
+module Sectors
+
+include("abstractcategory.jl")
+include("group_definitions.jl")
+include("category_definitions.jl")
+
+include("namedtuple_operations.jl")
+include("category_product.jl")
+
+end

NDTensors/src/lib/Sectors/src/abstractcategory.jl

@@ -0,0 +1,40 @@
+abstract type AbstractCategory end
+
+label(c::AbstractCategory) = error("method `label` not defined for type $(typeof(c))")
+
+function dimension(c::AbstractCategory)
+  return error("method `dimension` not defined for type $(typeof(c))")
+end
+
+function label_fusion_rule(::Type{C}, l1, l2) where {C<:AbstractCategory}
+  return error("label_fusion_rule not defined for type $(C)")
+end
+
+function fusion_rule(c1::AbstractCategory, c2::AbstractCategory)
+  C = typeof(c1)
+  return [C(d) for d in label_fusion_rule(C, label(c1), label(c2))]
+end
+
+⊗(c1::AbstractCategory, c2::AbstractCategory) = fusion_rule(c1, c2)
+
+Base.:(*)(s1::AbstractCategory, s2::AbstractCategory) = ⊗(s1, s2)
+
+⊕(a::AbstractCategory, b::AbstractCategory) = [a, b]
+⊕(v::Vector{<:AbstractCategory}, b::AbstractCategory) = vcat(v, b)
+⊕(a::AbstractCategory, v::Vector{<:AbstractCategory}) = vcat(a, v)
+
+function Base.show(io::IO, q::Vector{<:AbstractCategory})
+  (length(q) <= 1) && print(io, "[")
+  symbol = ""
+  for l in q
+    print(io, symbol, l)
+    symbol = " ⊕ "
+  end
+  return (length(q) <= 1) && print(io, "]")
+end
+
+function trivial(::Type{C}) where {C<:AbstractCategory}
+  return error("method `trivial` not defined for type $C")
+end
+
+istrivial(C::AbstractCategory) = (C == trivial(typeof(C)))

NDTensors/src/lib/Sectors/src/category_definitions.jl

@@ -0,0 +1,85 @@
+import HalfIntegers: Half
+
+#
+# Quantum 'group' su2ₖ
+#
+
+struct su2{k} <: AbstractCategory
+  j::Half{Int}
+end
+
+label(s::su2) = s.j
+
+level(s::su2{k}) where {k} = k
+
+trivial(::Type{su2{k}}) where {k} = su2{k}(0)
+
+function label_fusion_rule(::Type{su2{k}}, j1, j2) where {k}
+  return abs(j1 - j2):min(k - j1 - j2, j1 + j2)
+end
+
+#
+# Fibonacci category
+#
+# (same fusion rules as subcategory {0,1} of su2{3})
+#
+
+struct Fib <: AbstractCategory
+  l::Int
+end
+
+function Fib(s::AbstractString)
+  if s == "1"
+    return Fib(0)
+  elseif s == "τ"
+    return Fib(1)
+  end
+  return error("Unrecognized input \"$s\" to Fib constructor")
+end
+
+label(f::Fib) = f.l
+
+trivial(::Type{Fib}) = Fib(0)
+
+dimension(f::Fib) = istrivial(f) ? 1 : ((1 + √5) / 2)
+
+# Fusion rules identical to su2₃
+label_fusion_rule(::Type{Fib}, l1, l2) = label_fusion_rule(su2{3}, l1, l2)
+
+label_to_str(f::Fib) = istrivial(f) ? "1" : "τ"
+
+function Base.show(io::IO, f::Fib)
+  return print(io, "Fib(", label_to_str(f), ")")
+end
+
+#
+# Ising category
+#
+# (same fusion rules as su2{2})
+#
+
+struct Ising <: AbstractCategory
+  l::Half{Int}
+end
+
+function Ising(s::AbstractString)
+  for (a, v) in enumerate(("1", "σ", "ψ"))
+    (v == s) && return Ising((a - 1)//2)
+  end
+  return error("Unrecognized input \"$s\" to Ising constructor")
+end
+
+label(i::Ising) = i.l
+
+trivial(::Type{Ising}) = Ising(0)
+
+dimension(i::Ising) = (label(i) == 1//2) ? √2 : 1
+
+# Fusion rules identical to su2₂
+label_fusion_rule(::Type{Ising}, l1, l2) = label_fusion_rule(su2{2}, l1, l2)
+
+label_to_str(i::Ising) = ("1", "σ", "ψ")[Int(2 * label(i) + 1)]
+
+function Base.show(io::IO, f::Ising)
+  return print(io, "Ising(", label_to_str(f), ")")
+end

NDTensors/src/lib/Sectors/src/category_product.jl

@@ -0,0 +1,111 @@
+
+struct CategoryProduct{Categories} <: AbstractCategory
+  cats::Categories
+  global _CategoryProduct(l) = new{typeof(l)}(l)
+end
+
+CategoryProduct(c::CategoryProduct) = _CategoryProduct(categories(c))
+
+categories(s::CategoryProduct) = s.cats
+
+Base.isempty(S::CategoryProduct) = isempty(categories(S))
+Base.length(S::CategoryProduct) = length(categories(S))
+Base.getindex(S::CategoryProduct, args...) = getindex(categories(S), args...)
+
+function fusion_rule(s1::CategoryProduct, s2::CategoryProduct)
+  return [
+    CategoryProduct(l) for l in categories_fusion_rule(categories(s1), categories(s2))
+  ]
+end
+
+function Base.:(==)(A::CategoryProduct, B::CategoryProduct)
+  return categories_equal(categories(A), categories(B))
+end
+
+function Base.show(io::IO, s::CategoryProduct)
+  (length(s) < 2) && print(io, "sector")
+  print(io, "(")
+  symbol = ""
+  for p in pairs(categories(s))
+    print(io, symbol)
+    category_show(io, p[1], p[2])
+    symbol = " × "
+  end
+  return print(io, ")")
+end
+
+category_show(io::IO, k, v) = print(io, v)
+
+category_show(io::IO, k::Symbol, v) = print(io, "($k=$v,)")
+
+×(c1::AbstractCategory, c2::AbstractCategory) = ×(CategoryProduct(c1), CategoryProduct(c2))
+function ×(p1::CategoryProduct, p2::CategoryProduct)
+  return CategoryProduct(categories_product(categories(p1), categories(p2)))
+end
+
+categories_product(l1::NamedTuple, l2::NamedTuple) = union_keys(l1, l2)
+
+categories_product(l1::Tuple, l2::Tuple) = (l1..., l2...)
+
+×(nt1::NamedTuple, nt2::NamedTuple) = ×(CategoryProduct(nt1), CategoryProduct(nt2))
+×(c1::NamedTuple, c2::AbstractCategory) = ×(CategoryProduct(c1), CategoryProduct(c2))
+×(c1::AbstractCategory, c2::NamedTuple) = ×(CategoryProduct(c1), CategoryProduct(c2))
+
+#
+# Dictionary-like implementation
+#
+
+function CategoryProduct(nt::NamedTuple)
+  categories = sort_keys(nt)
+  return _CategoryProduct(categories)
+end
+
+CategoryProduct(; kws...) = CategoryProduct((; kws...))
+
+function CategoryProduct(pairs::Pair...)
+  keys = ntuple(n -> Symbol(pairs[n][1]), length(pairs))
+  vals = ntuple(n -> pairs[n][2], length(pairs))
+  return CategoryProduct(NamedTuple{keys}(vals))
+end
+
+function categories_fusion_rule(A::NamedTuple, B::NamedTuple)
+  qs = [A]
+  for (la, lb) in zip(pairs(intersect_keys(A, B)), pairs(intersect_keys(B, A)))
+    @assert la[1] == lb[1]
+    fused_vals = ⊗(la[2], lb[2])
+    qs = [union_keys((; la[1] => v), q) for v in fused_vals for q in qs]
+  end
+  # Include sectors of B not in A
+  qs = [union_keys(q, B) for q in qs]
+  return qs
+end
+
+function categories_equal(A::NamedTuple, B::NamedTuple)
+  common_categories = zip(pairs(intersect_keys(A, B)), pairs(intersect_keys(B, A)))
+  common_categories_match = all(nl -> (nl[1] == nl[2]), common_categories)
+  unique_categories_zero = all(l -> istrivial(l), symdiff_keys(A, B))
+  return common_categories_match && unique_categories_zero
+end
+
+#
+# Ordered implementation
+#
+
+CategoryProduct(t::Tuple) = _CategoryProduct(t)
+CategoryProduct(cats::AbstractCategory...) = CategoryProduct((cats...,))
+
+categories_equal(o1::Tuple, o2::Tuple) = (o1 == o2)
+
+function categories_fusion_rule(o1::Tuple, o2::Tuple)
+  N = length(o1)
+  length(o2) == N ||
+    throw(DimensionMismatch("Ordered CategoryProduct must have same size in ⊗"))
+  os = [o1]
+  replace(o, n, val) = ntuple(m -> (m == n) ? val : o[m], length(o))
+  for n in 1:N
+    os = [replace(o, n, f) for f in ⊗(o1[n], o2[n]) for o in os]
+  end
+  return os
+end
+
+sector(args...; kws...) = CategoryProduct(args...; kws...)