[lkdvos]

Braid and transpose a FusionTreeBlock of a ProductSector through its parts: the fusion trees of a product are pairs of fusion trees of two parts of the product, so the recoupling matrix of a block is the Kronecker product of the recoupling matrices of the two parts.

• A product with both UniqueFusion and other sectors is split into its other sectors and its UniqueFusion sectors. The latter have a single tree pair, such that a block is a block of the other sectors with a fixed tree pair and a scalar. All blocks that only differ in their abelian charges therefore share the cached manipulations of the non-abelian sectors: for fℤ₂ × U(1) × SU(2) × U(1) below, 26004 blocks map onto 142 SU(2) blocks.
• Any other product (e.g. several non-abelian sectors) is split into its first and its remaining sectors, recursively.
• UniqueFusion parts are transformed tree pair by tree pair directly in product_transform, so the generic fsbraid/fstranspose are left untouched.
• The canonical order of a product block, i.e. of FusionTreeBlock{I}(uncoupled, isdual) and therefore of fusionblocks, is now the Kronecker order of the two parts, such that the manipulations through the parts produce canonical blocks without reordering. This only affects the order of the trees within a block: the order of fusiontrees(::HomSpace), and thus the data layout of tensors, is unchanged.
• fusionblocks also benefits: the abelian part of a block is constructed first, which skips the non-abelian construction for the (many) combinations of uncoupled sectors with mismatching abelian charges.
• Named product sectors (NamedSector with at least two sectors) are handled in the same way, by deconstructing and reconstructing them as the corresponding ProductSector, so they share the caches of the unnamed sectors.

The code lives in src/fusiontrees/producttrees.jl: deconstructsector/reconstructsector and deconstructtree/reconstructtree split a sector or fusion tree into its two parts and join them back, productblock builds a block from the blocks of the parts, and product_transform combines their recoupling matrices. braid and transpose dispatch on ProductSectorBlock.

Cold construction of the recoupling data for one permutation of a tensor from a PEPS CTMRG step (the reproducers of #580 and #581), main → this PR, minimum of 5:

	main	this PR
fℤ₂ × U(1) × SU(3), 1 thread	6.2–6.7 s	3.7 s (1.7×)
fℤ₂ × U(1) × SU(3), 8 threads	1.7–2.0 s	1.0 s (1.7×)
fℤ₂ × U(1) × SU(2) × U(1), 1 thread	4.6–5.0 s	0.85 s (5.5×)
fℤ₂ × U(1) × SU(2) × U(1), 8 threads	2.6–3.1 s	0.75 s (3.8×)

Ranges are over two interleaved rounds; Julia 1.13.2, shared 32-core machine. For SU(2), fusionblocks (~0.4 s, single-threaded) is now the largest part of the remaining time, which is what #581 addresses; for SU(3), the remaining time is mostly braiding genuinely distinct SU(3) blocks, which is what #580 addresses. Both should compose with this PR.

Tested: the existing fusion tree tests (unitarity, reversibility, fusion tensor compatibility, planar traces) now exercise the manipulations through the parts for all product sectors in the test list. In that list, FermionParity ⊠ U1Irrep ⊠ SU2Irrep is replaced by a named Hubbard-like @NamedSector{f::FermionParity, s::SU2Irrep, c::U1Irrep} (the unnamed one stays in the fast list), which also exercises a non-abelian sector in between abelian ones. Additionally, for every block of a 3 ← 3 space of every product sector with a simple unit (including A4Irrep ⊠ A4Irrep with multiplicities, FermionParity ⊠ SU2Irrep ⊠ SU2Irrep and the anyonic products), and for the blocks of both reproducers, the results agree with the direct fsbraid/fstranspose to ~1e-15.

🤖 Generated with Claude Code

[leburgel]

This looks great, but I think there might be some additional gain for the specific case where there's only one symmetry factor which doesn't actually have UniqueFusion. This just happens to be the case I'm most interested in so maybe this would make things too artificial, but I think it could be worth it:

For products with a single factor without UniqueFusion (e.g. fℤ₂ × U(1) × SU(2) × U(1)), most of the time currently goes into the recursion around the braids rather than the braids themselves: splitting off the first factor gives three levels here, each copying all trees and building a block, and the result is then checked against the canonical block. All of this can be avoided by choosing the first factor without UniqueFusion as the head instead of the first factor, so that with a single such factor the tail is a product of UniqueFusion factors with a single tree per block. That needs three changes that build on each other:

1. Choose the head. Take the first factor without UniqueFusion. The factors before it have a single vertex label, so the head is still the slowest label of the product, as the splitting of the vertex labels in headtail assumes. With several such factors, the UniqueFusion factors then no longer get a level of their own either:

# blocks of products with `UniqueFusion` have a single tree and don't reach `product_transform` (see 3.)
@generated headfactor(::Type{ProductSector{T}}) where {T} =
    something(findfirst(F -> !(FusionStyle(F) isa UniqueFusion), collect(fieldtypes(T))), 1)
headsector(s::ProductSector) = s.sectors[headfactor(typeof(s))]
function tailsector(s::ProductSector)
    tail = _deletefactor(s.sectors, Val(headfactor(typeof(s))))
    return length(tail) == 1 ? only(tail) : ProductSector(tail)
end
@generated _deletefactor(t::Tuple, ::Val{k}) where {k} =
    Expr(:tuple, (:(t[$i]) for i in 1:fieldcount(t) if i != k)...)

2. Join head and tail at that position. The tree ⊠ puts the head first, so product_transform needs its inverse for the product type I instead (⊠ on trees isn't used elsewhere):

joinsectors(::Type{I}, h, t) where {I} =
    I(_insertfactor(t isa ProductSector ? t.sectors : (t,), h, Val(headfactor(I))))
@generated _insertfactor(t::Tuple, h, ::Val{k}) where {k} =
    Expr(:tuple, (i < k ? :(t[$i]) : i == k ? :h : :(t[$(i - 1)]) for i in 1:(fieldcount(t) + 1))...)

function jointrees(::Type{I}, h::FusionTree, t::FusionTree) where {I}
    vertices = ... # as in `⊠`, see 3.
    join(a, b) = joinsectors(I, a, b)
    return FusionTree(
        map(join, h.uncoupled, t.uncoupled), join(h.coupled, t.coupled), h.isdual,
        map(join, h.innerlines, t.innerlines), vertices
    )
end

with jointrees(I, h₁, t₁), jointrees(I, h₂, t₂) in place of h₁ ⊠ t₁, h₂ ⊠ t₂ in product_transform.

3. Use that the tail has unique fusion. With a single factor without UniqueFusion, the tail has UniqueFusion, which removes the remaining overhead in three places:

• A block of the tail has a single tree, so it doesn't need to recurse either: in braid(::ProductSectorBlock, ...) and transpose(::ProductSectorBlock, ...),

  (isempty(src) || FusionStyle(sectortype(src)) isa UniqueFusion) && return fsbraid((src, p, levels)) # resp. fstranspose

• The vertex labels of the tail are all 1, so the head keeps the labels of the product as they are, without vertexmultiplicities (and the same the other way around, in headtail and jointrees):

  # in `headtail`
  if FusionStyle(headsector(f.coupled)) isa UniqueFusion
      vertices_h = map(Returns(1), f.vertices)
  elseif FusionStyle(tailsector(f.coupled)) isa UniqueFusion
      vertices_t = map(Returns(1), f.vertices)
  else # as now
  ...
  # in `jointrees`
  if FusionStyle(I) isa GenericFusion && length(h) > 1
      vertices = FusionStyle(h) isa UniqueFusion ? t.vertices :
          FusionStyle(t) isa UniqueFusion ? h.vertices :
          map((μ, ν, n) -> ν + (μ - 1) * n, h.vertices, t.vertices, vertexmultiplicities(t))
  end

• With a single tree on one side, the product trees come out in the canonical order of the other side's block, which transform returns in canonical order, so the canonical block doesn't need to be rebuilt; in product_transform, before canonical = ...:

  (FusionStyle(first(first(src_h))) isa UniqueFusion || FusionStyle(first(first(src_t))) isa UniqueFusion) &&
      return dst => U

For products with several factors without UniqueFusion, the general head/tail recursion remains, now only over those factors. With these changes the braid transformers are identical to the current ones for fℤ₂ × U(1) × SU(3) and fℤ₂ × U(1) × SU(2) × U(1).

[leburgel]

One more thing for products with a single factor without UniqueFusion: the tail is now a product of UniqueFusion factors, whose blocks have a single tree, but those still go through product_transform factor by factor. They can go to fsbraid and fstranspose directly, which handle single trees already:

function braid(src::ProductSectorBlock, p::Index2Tuple, levels::Index2Tuple)
    # blocks of products without multiplicities have a single tree
    (isempty(src) || FusionStyle(sectortype(src)) isa UniqueFusion) && return fsbraid((src, p, levels))
    return product_transform(sectorscalartype(sectortype(src)), src) do x
        return braid(x, p, levels)
    end
end

and in the same way with fstranspose((src, p)) in transpose(src::ProductSectorBlock, p).

[lkdvos]

Thanks, these were very helpful! They are all in now, in a slightly different form:

1. Head choice. Rather than taking the first non-UniqueFusion factor as the head, a mixed product is now split into all its non-UniqueFusion factors (head) and all its UniqueFusion factors (tail), so the abelian factors never get a level of their own, also with several non-abelian factors. Only a product without UniqueFusion factors still recurses factor by factor.
2. Joining. joinsectors(I, h, t) interleaves the factors back into their original positions, and jointree(I, h, t) replaces the tree ⊠.
3. Unique tail.
   • With a UniqueFusion tail, product_transform transforms the single tail tree pair once, transforms the head block, and joins them, with U = coeff * U_head.
   • The vertex labels are not split or recombined in that case (no vertexmultiplicities).
   • The canonical order of a product block is now the Kronecker order itself, so there is no reordering into a canonical block at all anymore.
   • Blocks of all-UniqueFusion products are handled tree pair by tree pair inside product_transform instead of through fsbraid/fstranspose on the block. Doing it on the block tripped JET (fusiontrees of a tree pair) and fails for N = 0 blocks with several units.

With this, fℤ₂ × U(1) × SU(2) × U(1) goes from ~4.8 s to ~0.85 s cold (single-threaded); see the table in the description.

[lkdvos]

Alright, I think this one is ready for review. The main subtle change is that within a single block, the entries are now sorted such that the transformation matrix really is the kronecker product of the transformation matrix of the factors, which shouldn't change behavior at all but makes things both simpler and faster for this PR. This does not affect the data layout so in principle there is no issue with breaking changes etc.

[codecov]

Codecov Report

❌ Patch coverage is 93.69369% with 7 lines in your changes missing coverage. Please review.

Files with missing lines	Patch %	Lines
src/fusiontrees/producttrees.jl	94.00%	6 Missing ⚠️
ext/TensorKitEnzymeExt/utility.jl	0.00%	1 Missing ⚠️

Files with missing lines	Coverage Δ
src/fusiontrees/braiding_manipulations.jl	95.61% <100.00%> (+0.07%)	⬆️
src/fusiontrees/duality_manipulations.jl	97.15% <100.00%> (+0.02%)	⬆️
src/fusiontrees/fusiontrees.jl	96.15% <100.00%> (+0.04%)	⬆️
ext/TensorKitEnzymeExt/utility.jl	0.00% <0.00%> (ø)
src/fusiontrees/producttrees.jl	94.00% <94.00%> (ø)

... and 2 files with indirect coverage changes

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[Jutho]

Do we need this dichotomy ? What happens is otherfactors always returns a ProductSector? A ProductSector{Tuple{T}} with a single sector is not a ProductSectorType so it wouldn't be affected by anything in this file, no?

[lkdvos]

It's mostly an efficiency thing again: for example if we would have SU2Irrep boxtimes SU2Irrep, this would then be decomposed into SU2Irrep and ProductSector{Tuple{SU2Irrep}}, which then duplicates the cached entries we have for them. (also true for compile time)
In principle we can drop this, but I do think it is actually cleaner to have proper "factorization" that drops the ProductSector wrapper.

[Jutho]

Ok, I didn't manage to read through and grasp the complete logic yet. I find the naming somewhat confusing.

Firstly, this introduces the new term "factors" which is not very precise and clear. Secondly, there are quite a few names that sound very generic, such as unpartition, interleave, headtail, that might cause issues later since they are just in the general TensorKit namespace. If we need to keep them, I would definitely mark these as private methods by starting them with an underscore. I also find the naming a bit inconsistent. unpartition does at the ProductSector level what interleave does at the tuple level.

Then, firstfactor and otherfactors do exactly what first and tail do on tuples, but then at the level of sectors. So I would call these firstsector and tailsector.

The partitionfactor business is also very vague. Can we use something like deconstructsector, or, given that I assume the point is to split the sector into two parts, namely a non-unique and a unique fusion part, maybe bisectsector? Variable names g and a could perhaps become nu and u for non-unique and unique, or g and u for generic and unique (I don't think the fact that we often use u for unit is confusing here). The same naming convention but capitalized can be used for the corresponding types, which is currently also consistent.

For headsector and (current) tailsector, I am still thinking of what proper names could be. Checking whether I understand this correctly: in first instance, headsector will return the non-unique part, and only upon recursion, once the product sector is no longer mixed, it will actually return an individual (non-product) sector? Same for tailsector? So that is not completely in line with what to expect from head and tail. It is even not clear to me yet if this recursion even happens, or if you just want to split into the two parts with generic fusion and unique fusion, in all the headtail functionality on trees?

[lkdvos]

I'll give it another go. I agree the naming is a bit horrible, but I also have a hard time to come up with something better. To explain a bit more what is happening:

we are indeed factorizing a productsector into first all non-uniquefusion factors (which is where the recursion happens, these are kept separate), and then all uniquefusion factors (no recursion, these are all done in one go). Thus, this is not necessarily a split into two parts. There's also some subtle trickery that makes use of Trivial components because deligneproduct(x, Trivial() == x, which is how some of the additional factors end up being dropped.

Naming-wise, this is mostly because I originally factorized everything, which is when head and tail actually made sense, but I do agree now this is just misleading because that's not really what is happening anyways. Let me play around with some alternatives and push a new version in a bit.

[Jutho]

Maybe deconstructsector and reconstructsector can be useful terminology then.