From 6c14e697afd57fcc9b553ed6f6ce9e7da0af1f5c Mon Sep 17 00:00:00 2001 From: leburgel Date: Wed, 14 Oct 2026 12:59:06 +0200 Subject: [PATCH 1/7] Construct `FusionTreeBlock`s from the fusion product of the shorter side For sectors with a simple unit, a superset of the coupled sectors suffices, as the fusion tree iterators yield nothing for the others: the unit if one side is empty, and otherwise the fusion product of the shorter side only, instead of the intersection of the fusion products of both sides. The trees and their order are unchanged. Co-Authored-By: Claude Opus 5.5 --- src/fusiontrees/fusiontrees.jl | 20 +++++++++++++++----- 1 file changed, 15 insertions(+), 5 deletions(-) diff --git a/src/fusiontrees/fusiontrees.jl b/src/fusiontrees/fusiontrees.jl index 7b97223ec..7b2f289db 100644 --- a/src/fusiontrees/fusiontrees.jl +++ b/src/fusiontrees/fusiontrees.jl @@ -140,7 +140,15 @@ function FusionTreeBlock{I}( end # note: cs is unique due to ⊗ returning unique fusion results - cs = if N₁ == 0 + cs = if UnitStyle(I) isa SimpleUnit + # a superset of the coupled sectors suffices, as there are no fusion trees for the others: + # this avoids the fusion product of the longer side + if N₁ == 0 || N₂ == 0 + [unit(I)] + else + sort!(vec(collect(⊗(uncoupled[N₁ <= N₂ ? 1 : 2]...)))) # `vec`: products of `ProductSector`s have a shape + end + elseif N₁ == 0 sort!(collect(filter(isunit, ⊗(uncoupled[2]...)))) elseif N₂ == 0 sort!(collect(filter(isunit, ⊗(uncoupled[1]...)))) @@ -148,10 +156,12 @@ function FusionTreeBlock{I}( sort!(collect(intersect(⊗(uncoupled[1]...), ⊗(uncoupled[2]...)))) end - for c in cs, - f₂ in fusiontrees(uncoupled[2], c, isdual[2]), - f₁ in fusiontrees(uncoupled[1], c, isdual[1]) - push!(trees, (f₁, f₂)) + for c in cs + isempty(fusiontrees(uncoupled[1], c, isdual[1])) && continue + for f₂ in fusiontrees(uncoupled[2], c, isdual[2]), + f₁ in fusiontrees(uncoupled[1], c, isdual[1]) + push!(trees, (f₁, f₂)) + end end return FusionTreeBlock(trees) end From 2130dc3dad701113f76583818d2f6a2a0a5a5fc1 Mon Sep 17 00:00:00 2001 From: leburgel Date: Wed, 14 Oct 2026 12:59:06 +0200 Subject: [PATCH 2/7] Group the cached fusion trees in `fusionblocks` Instead of constructing a block for every combination of codomain and domain sectors, `fusionblocks` groups the cached `fusiontrees(W)` by their uncoupled sectors, in the same order as before. Co-Authored-By: Claude Opus 5.5 --- src/spaces/homspace.jl | 17 +++++++++++------ 1 file changed, 11 insertions(+), 6 deletions(-) diff --git a/src/spaces/homspace.jl b/src/spaces/homspace.jl index b41bd45be..45631db92 100644 --- a/src/spaces/homspace.jl +++ b/src/spaces/homspace.jl @@ -210,13 +210,18 @@ grouped by their uncoupled charges. function fusionblocks(W::HomSpace) I = sectortype(W) N₁, N₂ = numout(W), numin(W) - isdual_src = (map(isdual, codomain(W)), map(isdual, domain(W))) - fblocks = Vector{FusionTreeBlock{I, N₁, N₂, fusiontreetype(I, N₁, N₂)}}() - for dom_uncoupled_src in sectors(domain(W)), cod_uncoupled_src in sectors(codomain(W)) - fs_src = FusionTreeBlock{I}((cod_uncoupled_src, dom_uncoupled_src), isdual_src) - isempty(fs_src) || push!(fblocks, fs_src) + F = fusiontreetype(I, N₁, N₂) + # group the (cached) fusion trees of `W` by their uncoupled sectors + groups = Dict{Tuple{NTuple{N₁, I}, NTuple{N₂, I}}, Vector{F}}() + for (f₁, f₂) in fusiontrees(W) + push!(get!(Vector{F}, groups, (f₁.uncoupled, f₂.uncoupled)), (f₁, f₂)) end - return fblocks + # in the order of `sectors(domain(W))` (outer) and `sectors(codomain(W))` (inner), each of + # which runs fastest over its first space + position(V) = Dict{I, Int}(s => i for (i, s) in enumerate(sectors(V))) + codpos, dompos = map(position, codomain(W)), map(position, domain(W)) + rank((u₁, u₂)) = (reverse(map(getindex, dompos, u₂)), reverse(map(getindex, codpos, u₁))) + return map(u -> FusionTreeBlock(groups[u]), sort!(collect(keys(groups)); by = rank)) end function diagonalblockstructure(W::HomSpace) From 5d26ce33307a3d2de47d71da3ebc2ec0cd2d1a39 Mon Sep 17 00:00:00 2001 From: leburgel Date: Wed, 14 Oct 2026 14:47:29 +0200 Subject: [PATCH 3/7] Test the blocks of a HomSpace; changelog Co-Authored-By: Claude Opus 5.5 --- docs/src/Changelog.md | 1 + test/symmetries/doubletree.jl | 15 +++++++++++++++ 2 files changed, 16 insertions(+) diff --git a/docs/src/Changelog.md b/docs/src/Changelog.md index 2a853124a..b91e0e68b 100644 --- a/docs/src/Changelog.md +++ b/docs/src/Changelog.md @@ -42,6 +42,7 @@ When releasing a new version, move the "Unreleased" changes to a new version sec ### Performance - Recoupling matrices are stored such that they can leverage BLAS, also for mixed complex tensor with real recoupling coefficients. +- Constructing a `FusionTreeBlock` only takes the fusion product of the shorter side, and `fusionblocks` groups the cached fusion trees of the space instead of constructing a block for every combination of sectors (internal) ## [0.17.2](https://github.com/QuantumKitHub/TensorKit.jl/compare/v0.17.1...v0.17.2) - 2026-09-20 diff --git a/test/symmetries/doubletree.jl b/test/symmetries/doubletree.jl index 4fbe43e90..908e0497a 100644 --- a/test/symmetries/doubletree.jl +++ b/test/symmetries/doubletree.jl @@ -235,6 +235,21 @@ using TensorKitSectors end end + UnitStyle(I) isa SimpleUnit && @testset "Double fusion tree: blocks of a HomSpace" begin + V = Vect[I](Dict(s => 1 for s in (out..., incoming))) + W = V ⊗ V' ← V ⊗ V' # never empty: the unit couples both sides + trees = collect(fusiontrees(W)) + blocks = TK.fusionblocks(W) + # every block holds the trees of `W` with its uncoupled sectors, in the same order + for block in blocks + @test fusiontrees(block) == filter(f -> (f[1].uncoupled, f[2].uncoupled) == block.uncoupled, trees) + end + @test sum(length, blocks) == length(trees) + # the blocks are ordered as `sectors(domain(W))` (outer) and `sectors(codomain(W))` (inner) + uncoupled = [(u₁, u₂) for u₂ in sectors(domain(W)) for u₁ in sectors(codomain(W))] + @test map(block -> block.uncoupled, blocks) == filter(u -> any(f -> (f[1].uncoupled, f[2].uncoupled) == u, trees), uncoupled) + end + @testset "Double fusion tree: planar trace" begin if FusionStyle(I) isa UniqueFusion f1, f1 = src From c643ae2151e8cfec693c960e561b682457293f75 Mon Sep 17 00:00:00 2001 From: leburgel Date: Thu, 15 Oct 2026 13:22:59 +0200 Subject: [PATCH 4/7] Compute the rank of each block once when sorting in `fusionblocks` `sort!(...; by = rank)` evaluated `rank`, and with it the `Dict` lookups of (product) sectors, at every comparison. Co-Authored-By: Claude Opus 5.5 --- src/spaces/homspace.jl | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/spaces/homspace.jl b/src/spaces/homspace.jl index 45631db92..dc5bd589b 100644 --- a/src/spaces/homspace.jl +++ b/src/spaces/homspace.jl @@ -221,7 +221,8 @@ function fusionblocks(W::HomSpace) position(V) = Dict{I, Int}(s => i for (i, s) in enumerate(sectors(V))) codpos, dompos = map(position, codomain(W)), map(position, domain(W)) rank((u₁, u₂)) = (reverse(map(getindex, dompos, u₂)), reverse(map(getindex, codpos, u₁))) - return map(u -> FusionTreeBlock(groups[u]), sort!(collect(keys(groups)); by = rank)) + us = collect(keys(groups)) + return map(i -> FusionTreeBlock(groups[us[i]]), sortperm(map(rank, us))) # each rank computed once end function diagonalblockstructure(W::HomSpace) From 5837d51c5f333aa2c594fb33d059f3d5040a3ebd Mon Sep 17 00:00:00 2001 From: leburgel Date: Thu, 15 Oct 2026 13:25:51 +0200 Subject: [PATCH 5/7] Use the fusion product of the shorter side for every UnitStyle; check only the longer side for trees Co-Authored-By: Claude Opus 5.5 --- src/fusiontrees/fusiontrees.jl | 23 +++++++++-------------- 1 file changed, 9 insertions(+), 14 deletions(-) diff --git a/src/fusiontrees/fusiontrees.jl b/src/fusiontrees/fusiontrees.jl index 7b2f289db..8684f8371 100644 --- a/src/fusiontrees/fusiontrees.jl +++ b/src/fusiontrees/fusiontrees.jl @@ -140,24 +140,19 @@ function FusionTreeBlock{I}( end # note: cs is unique due to ⊗ returning unique fusion results - cs = if UnitStyle(I) isa SimpleUnit - # a superset of the coupled sectors suffices, as there are no fusion trees for the others: - # this avoids the fusion product of the longer side - if N₁ == 0 || N₂ == 0 - [unit(I)] - else - sort!(vec(collect(⊗(uncoupled[N₁ <= N₂ ? 1 : 2]...)))) # `vec`: products of `ProductSector`s have a shape - end - elseif N₁ == 0 - sort!(collect(filter(isunit, ⊗(uncoupled[2]...)))) - elseif N₂ == 0 - sort!(collect(filter(isunit, ⊗(uncoupled[1]...)))) + # a superset of the coupled sectors suffices, as there are no fusion trees for the others: + # the units, or the fusion product of the shorter side, which avoids that of the longer side + # (`vec`: the units and products of `ProductSector`s have a shape) + cs = if N₁ == 0 || N₂ == 0 + sort!(vec(collect(allunits(I)))) else - sort!(collect(intersect(⊗(uncoupled[1]...), ⊗(uncoupled[2]...)))) + sort!(vec(collect(⊗(uncoupled[N₁ <= N₂ ? 1 : 2]...)))) end for c in cs - isempty(fusiontrees(uncoupled[1], c, isdual[1])) && continue + # the shorter side has fusion trees for every candidate, the longer side not necessarily + longer = N₁ <= N₂ ? fusiontrees(uncoupled[2], c, isdual[2]) : fusiontrees(uncoupled[1], c, isdual[1]) + isempty(longer) && continue for f₂ in fusiontrees(uncoupled[2], c, isdual[2]), f₁ in fusiontrees(uncoupled[1], c, isdual[1]) push!(trees, (f₁, f₂)) From 2cb2bf70e319d9bc0ef8a91d59f878d30712497a Mon Sep 17 00:00:00 2001 From: leburgel Date: Thu, 15 Oct 2026 20:37:14 +0200 Subject: [PATCH 6/7] Use the unit of the other side as the coupled sector of a block with no legs on one side Co-Authored-By: Claude Opus 5.5 --- src/fusiontrees/fusiontrees.jl | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/src/fusiontrees/fusiontrees.jl b/src/fusiontrees/fusiontrees.jl index 8684f8371..45f2c80a4 100644 --- a/src/fusiontrees/fusiontrees.jl +++ b/src/fusiontrees/fusiontrees.jl @@ -141,10 +141,12 @@ function FusionTreeBlock{I}( # note: cs is unique due to ⊗ returning unique fusion results # a superset of the coupled sectors suffices, as there are no fusion trees for the others: - # the units, or the fusion product of the shorter side, which avoids that of the longer side - # (`vec`: the units and products of `ProductSector`s have a shape) - cs = if N₁ == 0 || N₂ == 0 - sort!(vec(collect(allunits(I)))) + # the unit of the other side, or the fusion product of the shorter side, which avoids that of + # the longer side (`vec`: products of `ProductSector`s have a shape) + cs = if N₁ == 0 + [leftunit(first(uncoupled[2]))] + elseif N₂ == 0 + [leftunit(first(uncoupled[1]))] else sort!(vec(collect(⊗(uncoupled[N₁ <= N₂ ? 1 : 2]...)))) end From ed344a52781415a3ccc7848ac970be4d5ad079e6 Mon Sep 17 00:00:00 2001 From: leburgel Date: Thu, 15 Oct 2026 22:21:32 +0200 Subject: [PATCH 7/7] Only try coupled sectors whose units match the ends of the other side For several units, a coupled sector of a block must be compatible with the units at both ends of either side; otherwise the fusion tree iterator throws instead of being empty. Co-Authored-By: Claude Opus 5.5 --- src/fusiontrees/fusiontrees.jl | 12 ++++++++---- 1 file changed, 8 insertions(+), 4 deletions(-) diff --git a/src/fusiontrees/fusiontrees.jl b/src/fusiontrees/fusiontrees.jl index 45f2c80a4..093d9daca 100644 --- a/src/fusiontrees/fusiontrees.jl +++ b/src/fusiontrees/fusiontrees.jl @@ -143,16 +143,20 @@ function FusionTreeBlock{I}( # a superset of the coupled sectors suffices, as there are no fusion trees for the others: # the unit of the other side, or the fusion product of the shorter side, which avoids that of # the longer side (`vec`: products of `ProductSector`s have a shape) - cs = if N₁ == 0 - [leftunit(first(uncoupled[2]))] - elseif N₂ == 0 - [leftunit(first(uncoupled[1]))] + cs = if N₁ == 0 || N₂ == 0 + # the only unit that can couple: on the left of the first and on the right of the last sector + us = N₁ == 0 ? uncoupled[2] : uncoupled[1] + leftunit(first(us)) == rightunit(last(us)) ? [leftunit(first(us))] : I[] else sort!(vec(collect(⊗(uncoupled[N₁ <= N₂ ? 1 : 2]...)))) end for c in cs # the shorter side has fusion trees for every candidate, the longer side not necessarily + # (with several units, `c` must also have the units at the ends of the longer side) + us = N₁ <= N₂ ? uncoupled[2] : uncoupled[1] + UnitStyle(I) isa SimpleUnit || N₁ == 0 || N₂ == 0 || + (leftunit(c) == leftunit(first(us)) && rightunit(c) == rightunit(last(us))) || continue longer = N₁ <= N₂ ? fusiontrees(uncoupled[2], c, isdual[2]) : fusiontrees(uncoupled[1], c, isdual[1]) isempty(longer) && continue for f₂ in fusiontrees(uncoupled[2], c, isdual[2]),