diff --git a/Project.toml b/Project.toml index ee0408a..250bfd2 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "GradedUnitRanges" uuid = "e2de450a-8a67-46c7-b59c-01d5a3d041c5" authors = ["ITensor developers and contributors"] -version = "0.1.7" +version = "0.2.0" [deps] BlockArrays = "8e7c35d0-a365-5155-bbbb-fb81a777f24e" @@ -9,6 +9,7 @@ Compat = "34da2185-b29b-5c13-b0c7-acf172513d20" FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b" LabelledNumbers = "f856a3a6-4152-4ec4-b2a7-02c1a55d7993" SplitApplyCombine = "03a91e81-4c3e-53e1-a0a4-9c0c8f19dd66" +TensorProducts = "decf83d6-1968-43f4-96dc-fdb3fe15fc6d" [weakdeps] SymmetrySectors = "f8a8ad64-adbc-4fce-92f7-ffe2bb36a86e" @@ -23,4 +24,5 @@ FillArrays = "1.13.0" LabelledNumbers = "0.1.0" SplitApplyCombine = "1.2.3" SymmetrySectors = "0.1.4" +TensorProducts = "0.1.2" julia = "1.10" diff --git a/docs/Project.toml b/docs/Project.toml index 6e69b68..9f52bfe 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -6,6 +6,6 @@ Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306" [compat] Documenter = "1.8.1" -GradedUnitRanges = "0.1.6" +GradedUnitRanges = "0.2" LabelledNumbers = "0.1.0" Literate = "2.20.1" diff --git a/src/GradedUnitRanges.jl b/src/GradedUnitRanges.jl index e675bc4..90fb410 100644 --- a/src/GradedUnitRanges.jl +++ b/src/GradedUnitRanges.jl @@ -7,7 +7,6 @@ include("gradedunitrange.jl") include("dual.jl") include("labelledunitrangedual.jl") include("gradedunitrangedual.jl") -include("onetoone.jl") include("fusion.jl") end diff --git a/src/fusion.jl b/src/fusion.jl index 14ea931..2a18bf0 100644 --- a/src/fusion.jl +++ b/src/fusion.jl @@ -1,40 +1,10 @@ using BlockArrays: AbstractBlockedUnitRange, blocklengths using LabelledNumbers: LabelledInteger, label, labelled using SplitApplyCombine: groupcount - -# https://github.com/ITensor/ITensors.jl/blob/v0.3.57/NDTensors/src/lib/GradedAxes/src/tensor_product.jl -# https://en.wikipedia.org/wiki/Tensor_product -# https://github.com/KeitaNakamura/Tensorial.jl -function tensor_product( - a1::AbstractUnitRange, - a2::AbstractUnitRange, - a3::AbstractUnitRange, - a_rest::Vararg{AbstractUnitRange}, -) - return foldl(tensor_product, (a1, a2, a3, a_rest...)) -end +using TensorProducts: TensorProducts, OneToOne, tensor_product flip_dual(r::AbstractUnitRange) = r flip_dual(r::GradedUnitRangeDual) = flip(r) -function tensor_product(a1::AbstractUnitRange, a2::AbstractUnitRange) - return tensor_product(flip_dual(a1), flip_dual(a2)) -end - -function tensor_product(a1::Base.OneTo, a2::Base.OneTo) - return Base.OneTo(length(a1) * length(a2)) -end - -function tensor_product(::OneToOne, a2::AbstractUnitRange) - return a2 -end - -function tensor_product(a1::AbstractUnitRange, ::OneToOne) - return a1 -end - -function tensor_product(::OneToOne, ::OneToOne) - return OneToOne() -end function fuse_labels(x, y) return error( @@ -42,17 +12,19 @@ function fuse_labels(x, y) ) end -function fuse_blocklengths(x::Integer, y::Integer) - # return blocked unit range to keep non-abelian interface - return blockedrange([x * y]) -end - function fuse_blocklengths(x::LabelledInteger, y::LabelledInteger) # return blocked unit range to keep non-abelian interface return blockedrange([labelled(x * y, fuse_labels(label(x), label(y)))]) end -function tensor_product(a1::AbstractBlockedUnitRange, a2::AbstractBlockedUnitRange) +unmerged_tensor_product() = OneToOne() +unmerged_tensor_product(a) = a +unmerged_tensor_product(a1, a2) = tensor_product(a1, a2) +function unmerged_tensor_product(a1, a2, as...) + return unmerged_tensor_product(unmerged_tensor_product(a1, a2), as...) +end + +function unmerged_tensor_product(a1::AbstractGradedUnitRange, a2::AbstractGradedUnitRange) nested = map(Iterators.flatten((Iterators.product(blocks(a1), blocks(a2)),))) do it return mapreduce(length, fuse_blocklengths, it) end @@ -96,15 +68,11 @@ end blockmergesort(g::GradedUnitRangeDual) = flip(blockmergesort(flip(g))) blockmergesort(g::AbstractUnitRange) = g -# fusion_product produces a sorted, non-dual GradedUnitRange -function fusion_product(g1, g2) - return blockmergesort(tensor_product(g1, g2)) -end +# tensor_product produces a sorted, non-dual GradedUnitRange +TensorProducts.tensor_product(g::AbstractGradedUnitRange) = blockmergesort(flip_dual(g)) -fusion_product(g::AbstractUnitRange) = blockmergesort(g) -fusion_product(g::GradedUnitRangeDual) = fusion_product(flip(g)) - -# recursive fusion_product. Simpler than reduce + fix type stability issues with reduce -function fusion_product(g1, g2, g3...) - return fusion_product(fusion_product(g1, g2), g3...) +function TensorProducts.tensor_product( + g1::AbstractGradedUnitRange, g2::AbstractGradedUnitRange +) + return blockmergesort(unmerged_tensor_product(g1, g2)) end diff --git a/src/onetoone.jl b/src/onetoone.jl deleted file mode 100644 index 3f7b999..0000000 --- a/src/onetoone.jl +++ /dev/null @@ -1,8 +0,0 @@ -using BlockArrays: BlockArrays - -# Represents the range `1:1` or `Base.OneTo(1)`. -struct OneToOne{T} <: AbstractUnitRange{T} end -OneToOne() = OneToOne{Bool}() -Base.first(a::OneToOne) = one(eltype(a)) -Base.last(a::OneToOne) = one(eltype(a)) -BlockArrays.blockaxes(g::OneToOne) = (Block.(g),) # BlockArrays default crashes for OneToOne{Bool} diff --git a/test/Project.toml b/test/Project.toml index 4ba9a86..c1f51f5 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -1,19 +1,19 @@ [deps] Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595" BlockArrays = "8e7c35d0-a365-5155-bbbb-fb81a777f24e" -BlockSparseArrays = "2c9a651f-6452-4ace-a6ac-809f4280fbb4" GradedUnitRanges = "e2de450a-8a67-46c7-b59c-01d5a3d041c5" LabelledNumbers = "f856a3a6-4152-4ec4-b2a7-02c1a55d7993" SafeTestsets = "1bc83da4-3b8d-516f-aca4-4fe02f6d838f" Suppressor = "fd094767-a336-5f1f-9728-57cf17d0bbfb" +TensorProducts = "decf83d6-1968-43f4-96dc-fdb3fe15fc6d" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [compat] Aqua = "0.8.9" BlockArrays = "1.4.0" -BlockSparseArrays = "0.2.28, 0.3" -GradedUnitRanges = "0.1.6" +GradedUnitRanges = "0.2" LabelledNumbers = "0.1.0" SafeTestsets = "0.1" Suppressor = "0.2" +TensorProducts = "0.1.0" Test = "1.10" diff --git a/test/test_basics.jl b/test/test_basics.jl index cbd9bc4..bae6eb0 100644 --- a/test/test_basics.jl +++ b/test/test_basics.jl @@ -1,4 +1,3 @@ -@eval module $(gensym()) using BlockArrays: Block, BlockRange, @@ -13,34 +12,13 @@ using BlockArrays: combine_blockaxes, mortar using GradedUnitRanges: - GradedOneTo, - GradedUnitRange, - OneToOne, - blocklabels, - gradedrange, - sector_type, - space_isequal + GradedOneTo, GradedUnitRange, blocklabels, gradedrange, sector_type, space_isequal using LabelledNumbers: LabelledUnitRange, islabelled, label, labelled, labelled_isequal, unlabel using Test: @test, @test_broken, @testset -@testset "OneToOne" begin - a0 = OneToOne() - @test a0 isa OneToOne{Bool} - @test eltype(a0) == Bool - @test length(a0) == 1 - @test labelled_isequal(a0, a0) - @test a0[1] == true - @test a0[[1]] == [true] - - @test labelled_isequal(a0, 1:1) - @test labelled_isequal(1:1, a0) - @test !labelled_isequal(a0, 1:2) - @test !labelled_isequal(1:2, a0) -end - @testset "GradedUnitRanges basics" begin - a0 = OneToOne() + a0 = Base.OneTo(1) for a in ( blockedrange([labelled(2, "x"), labelled(3, "y")]), gradedrange([labelled(2, "x"), labelled(3, "y")]), @@ -260,4 +238,3 @@ end @test length(a) == 1 @test label(first(a)) == "x" end -end diff --git a/test/test_dual.jl b/test/test_dual.jl index 7feb159..2ac6777 100644 --- a/test/test_dual.jl +++ b/test/test_dual.jl @@ -1,4 +1,3 @@ -@eval module $(gensym()) using BlockArrays: Block, BlockedOneTo, @@ -13,13 +12,11 @@ using BlockArrays: findblock, mortar, combine_blockaxes -using BlockSparseArrays: BlockSparseArray using GradedUnitRanges: AbstractGradedUnitRange, GradedUnitRanges, GradedUnitRangeDual, LabelledUnitRangeDual, - OneToOne, blocklabels, blockmergesortperm, blocksortperm, @@ -36,6 +33,8 @@ using GradedUnitRanges: using LabelledNumbers: LabelledInteger, LabelledUnitRange, label, label_type, labelled, labelled_isequal, unlabel using Test: @test, @test_broken, @testset +using TensorProducts: OneToOne, tensor_product + struct U1 n::Int end @@ -306,16 +305,3 @@ end @test !isdual(dual(flip(a))) end end - -@testset "dag" begin - elt = ComplexF64 - r = gradedrange([U1(0) => 2, U1(1) => 3]) - a = BlockSparseArray{elt}(undef, r, dual(r)) - a[Block(1, 1)] = randn(elt, 2, 2) - a[Block(2, 2)] = randn(elt, 3, 3) - @test isdual.(axes(a)) == (false, true) - ad = dag(a) - @test Array(ad) == conj(Array(a)) - @test isdual.(axes(ad)) == (true, false) -end -end diff --git a/test/test_tensor_product.jl b/test/test_tensor_product.jl index 604f1d0..f58a447 100644 --- a/test/test_tensor_product.jl +++ b/test/test_tensor_product.jl @@ -1,4 +1,3 @@ -@eval module $(gensym()) using Test: @test, @testset using BlockArrays: blocklength, blocklengths @@ -6,16 +5,15 @@ using BlockArrays: blocklength, blocklengths using GradedUnitRanges: GradedUnitRanges, GradedOneTo, - OneToOne, + blocklabels, dual, - fusion_product, + unmerged_tensor_product, flip, gradedrange, space_isequal, - isdual, - tensor_product - + isdual using LabelledNumbers: labelled_isequal +using TensorProducts: OneToOne, tensor_product struct U1 n::Int @@ -23,66 +21,90 @@ end GradedUnitRanges.dual(c::U1) = U1(-c.n) Base.isless(c1::U1, c2::U1) = c1.n < c2.n GradedUnitRanges.fuse_labels(x::U1, y::U1) = U1(x.n + y.n) +a0 = gradedrange([U1(1) => 1, U1(2) => 3, U1(1) => 1]) -@testset "GradedUnitRanges.tensor_product" begin +@testset "unmerged_tensor_product" begin GradedUnitRanges.fuse_labels(x::String, y::String) = x * y - g0 = OneToOne() - @test labelled_isequal(g0, g0) - @test labelled_isequal(tensor_product(g0, g0), g0) + @test unmerged_tensor_product() isa OneToOne a = gradedrange(["x" => 2, "y" => 3]) - b = tensor_product(a, a) + @test labelled_isequal(unmerged_tensor_product(a), a) + + b = unmerged_tensor_product(a, a) @test b isa GradedOneTo @test length(b) == 25 @test blocklength(b) == 4 @test blocklengths(b) == [4, 6, 6, 9] @test labelled_isequal(b, gradedrange(["xx" => 4, "yx" => 6, "xy" => 6, "yy" => 9])) - c = tensor_product(a, a, a) + c = unmerged_tensor_product(a, a, a) @test c isa GradedOneTo @test length(c) == 125 @test blocklength(c) == 8 -end - -@testset "GradedUnitRanges.fusion_product" begin - g0 = OneToOne() - @test labelled_isequal(fusion_product(g0, g0), g0) + @test blocklabels(c) == ["xxx", "yxx", "xyx", "yyx", "xxy", "yxy", "xyy", "yyy"] - a = gradedrange([U1(1) => 1, U1(2) => 3, U1(1) => 1]) - - b = fusion_product(a) - @test labelled_isequal(b, gradedrange([U1(1) => 2, U1(2) => 3])) - - c = fusion_product(a, a) - @test labelled_isequal(c, gradedrange([U1(2) => 4, U1(3) => 12, U1(4) => 9])) - - d = fusion_product(a, a, a) @test labelled_isequal( - d, gradedrange([U1(3) => 8, U1(4) => 36, U1(5) => 54, U1(6) => 27]) + unmerged_tensor_product(a0, a0), + gradedrange([ + U1(2) => 1, + U1(3) => 3, + U1(2) => 1, + U1(3) => 3, + U1(4) => 9, + U1(3) => 3, + U1(2) => 1, + U1(3) => 3, + U1(2) => 1, + ]), ) end -@testset "dual and tensor_product" begin - a = gradedrange([U1(1) => 1, U1(2) => 3, U1(1) => 1]) - ad = dual(a) - - b = fusion_product(ad) - @test b isa GradedOneTo - @test !isdual(b) - @test space_isequal(b, gradedrange([U1(-2) => 3, U1(-1) => 2])) - - c = fusion_product(ad, ad) - @test c isa GradedOneTo - @test !isdual(c) - @test space_isequal(c, gradedrange([U1(-4) => 9, U1(-3) => 12, U1(-2) => 4])) - - d = fusion_product(ad, a) - @test !isdual(d) - @test space_isequal(d, gradedrange([U1(-1) => 6, U1(0) => 13, U1(1) => 6])) - - e = fusion_product(a, ad) - @test !isdual(d) - @test space_isequal(e, d) +@testset "symmetric tensor_product" begin + for a in (a0, a0[1:5]) + @test labelled_isequal(tensor_product(a), gradedrange([U1(1) => 2, U1(2) => 3])) + + @test labelled_isequal( + tensor_product(a, a), gradedrange([U1(2) => 4, U1(3) => 12, U1(4) => 9]) + ) + @test labelled_isequal( + tensor_product(a, OneToOne()), gradedrange([U1(1) => 2, U1(2) => 3]) + ) + @test labelled_isequal( + tensor_product(OneToOne(), a), gradedrange([U1(1) => 2, U1(2) => 3]) + ) + + d = tensor_product(a, a, a) + @test labelled_isequal( + d, gradedrange([U1(3) => 8, U1(4) => 36, U1(5) => 54, U1(6) => 27]) + ) + end end + +@testset "dual and tensor_product" begin + for a in (a0, a0[1:5]) + ad = dual(a) + + b = unmerged_tensor_product(ad) + @test isdual(b) + @test space_isequal(b, ad) + + b = tensor_product(ad) + @test b isa GradedOneTo + @test !isdual(b) + @test space_isequal(b, gradedrange([U1(-2) => 3, U1(-1) => 2])) + + c = tensor_product(ad, ad) + @test c isa GradedOneTo + @test !isdual(c) + @test space_isequal(c, gradedrange([U1(-4) => 9, U1(-3) => 12, U1(-2) => 4])) + + d = tensor_product(ad, a) + @test !isdual(d) + @test space_isequal(d, gradedrange([U1(-1) => 6, U1(0) => 13, U1(1) => 6])) + + e = tensor_product(a, ad) + @test !isdual(d) + @test space_isequal(e, d) + end end