diff --git a/src/tensors/abstracttensor.jl b/src/tensors/abstracttensor.jl index 48bfb0dd7..839b260bf 100644 --- a/src/tensors/abstracttensor.jl +++ b/src/tensors/abstracttensor.jl @@ -130,6 +130,7 @@ Determine an appropriate storage type for the combination of tensors `A` and `B` Optionally, a scalartype `T` for the destination can be supplied that might differ from the inputs. """ promote_storagetype +@inline promote_storagetype(A::AbstractTensorMap) = storagetype(A) @inline promote_storagetype(A::AbstractTensorMap, B::AbstractTensorMap, Cs::AbstractTensorMap...) = promote_storagetype(storagetype(A), storagetype(B), map(storagetype, Cs)...) @inline promote_storagetype(::Type{T}, A::AbstractTensorMap, B::AbstractTensorMap, Cs::AbstractTensorMap...) where {T <: Number} = diff --git a/src/tensors/braidingtensor.jl b/src/tensors/braidingtensor.jl index f08dd8181..d28b2e1df 100644 --- a/src/tensors/braidingtensor.jl +++ b/src/tensors/braidingtensor.jl @@ -199,121 +199,113 @@ function add_transform!( ) end -# VectorInterface -# --------------- -# TODO - -# TensorOperations -# ---------------- -# TODO: implement specialized methods - -function TO.tensoradd!( - C::AbstractTensorMap, - A::BraidingTensor, pA::Index2Tuple, conjA::Symbol, - α::Number, β::Number, backend = TO.DefaultBackend(), allocator = TO.DefaultAllocator() - ) - return TO.tensoradd!(C, TensorMap(A), pA, conjA, α, β, backend, allocator) -end - -# Planar operations -# ----------------- -# TODO: implement specialized methods - -function planaradd!( - C::AbstractTensorMap, - A::BraidingTensor, p::Index2Tuple, - α::Number, β::Number, backend, allocator - ) - return planaradd!(C, TensorMap(A), p, α, β, backend, allocator) -end - function planarcontract!( C::AbstractTensorMap, - A::BraidingTensor, - (oindA, cindA)::Index2Tuple, - B::AbstractTensorMap, - (cindB, oindB)::Index2Tuple, - (p1, p2)::Index2Tuple, + A::BraidingTensor, pA::Index2Tuple, + B::AbstractTensorMap, pB::Index2Tuple, + pAB::Index2Tuple, α::Number, β::Number, backend, allocator ) # special case only defined for contracting 2 indices - length(oindA) == length(cindA) == 2 || - return planarcontract!( - C, TensorMap(A), (oindA, cindA), B, (cindB, oindB), (p1, p2), - α, β, backend, allocator - ) + length.(pA) == (2, 2) || + return planarcontract!(C, TensorMap(A), pA, B, pB, pAB, α, β, backend, allocator) + + spacecheck_contract(C, A, pA, false, B, pB, false, pAB) codA, domA = codomainind(A), domainind(A) codB, domB = codomainind(B), domainind(B) oindA, cindA, oindB, cindB = reorder_indices( - codA, domA, codB, domB, oindA, cindA, oindB, cindB, p1, p2 + codA, domA, codB, domB, pA..., reverse(pB)..., pAB... ) - if space(B, cindB[1]) != space(A, cindA[1])' || - space(B, cindB[2]) != space(A, cindA[2])' - throw(SpaceMismatch("$(space(C)) ≠ permute($(space(A))[$oindA, $cindA] * $(space(B))[$cindB, $oindB], ($p1, $p2)")) - end - - if BraidingStyle(sectortype(B)) isa Bosonic + I = sectortype(C) + BraidingStyle(I) isa Bosonic && return add_permute!(C, B, (reverse(cindB), oindB), α, β, backend) - end - τ_levels = A.adjoint ? (1, 2, 2, 1) : (2, 1, 1, 2) - scale!(C, β) - - inv_braid = τ_levels[cindA[1]] > τ_levels[cindA[2]] - for (f₁, f₂) in fusiontrees(B) - local newtrees - for ((f₁′, f₂′), coeff′) in transpose((f₁, f₂), (cindB, oindB)) - for (f₁′′, coeff′′) in artin_braid(f₁′, 1; inv = inv_braid) - f12 = (f₁′′, f₂′) - coeff = coeff′ * coeff′′ - if @isdefined newtrees - newtrees[f12] = get(newtrees, f12, zero(coeff)) + coeff - else - newtrees = Dict(f12 => coeff) - end - end - end - for ((f₁′, f₂′), coeff) in newtrees - TO.tensoradd!( - C[f₁′, f₂′], B[f₁, f₂], (reverse(cindB), oindB), false, - α * coeff, One(), backend, allocator - ) - end + # Non-bosonic case: factor into a cyclic transpose (no crossings) + a single Artin braid + # that swaps the two contracted legs, producing the R-symbol that A encodes. Naively + # using a single `add_braid!` is wrong: it would resolve cyclic moves as crossings and + # pick up spurious R-symbol factors. + B_in_layout = (cindB == codB && oindB == domB) + if B_in_layout + B′ = B + else + B′ = TO.tensoralloc_add( + scalartype(B), B, (cindB, oindB), false, Val(true), allocator + ) + add_transpose!(B′, B, (cindB, oindB), One(), Zero(), backend) end + + levelsA = A.adjoint ? (1, 2, 2, 1) : (2, 1, 1, 2) + N = numind(B) + levels = ( + levelsA[cindA[1]], levelsA[cindA[2]], + ntuple(Returns(3), N - 2)..., + ) + + add_braid!( + C, B′, ((2, 1), ntuple(i -> i + 2, N - 2)), + levels, α, β, backend, + ) + + B_in_layout || TO.tensorfree!(B′, allocator) return C end function planarcontract!( C::AbstractTensorMap, - A::AbstractTensorMap, (oindA, cindA)::Index2Tuple, - B::BraidingTensor, (cindB, oindB)::Index2Tuple, - (p1, p2)::Index2Tuple, + A::AbstractTensorMap, pA::Index2Tuple, + B::BraidingTensor, pB::Index2Tuple, + pAB::Index2Tuple, α::Number, β::Number, backend, allocator ) - # special case only defined for contracting 2 indices - length(oindB) == length(cindB) == 2 || - return planarcontract!( - C, A, (oindA, cindA), TensorMap(B), (cindB, oindB), (p1, p2), - α, β, backend, allocator - ) + # special case only defined for contracting all 4 indices of B (2 contracted + 2 open) + length.(pB) == (2, 2) || + return planarcontract!(C, A, pA, TensorMap(B), pB, pAB, α, β, backend, allocator) + + spacecheck_contract(C, A, pA, false, B, pB, false, pAB) codA, domA = codomainind(A), domainind(A) codB, domB = codomainind(B), domainind(B) oindA, cindA, oindB, cindB = reorder_indices( - codA, domA, codB, domB, oindA, cindA, oindB, cindB, p1, p2 + codA, domA, codB, domB, pA..., reverse(pB)..., pAB... ) - if space(B, cindB[1]) != space(A, cindA[1])' || space(B, cindB[2]) != space(A, cindA[2])' - throw(SpaceMismatch("$(space(C)) ≠ permute($(space(A))[$oindA, $cindA] * $(space(B))[$cindB, $oindB], ($p1, $p2)")) + I = sectortype(C) + BraidingStyle(I) isa Bosonic && + return add_permute!(C, A, (oindA, reverse(cindA)), α, β, backend) + + # Non-bosonic case: cyclic transpose A → (oindA, cindA) (no crossings), then a single + # Artin braid swaps A′'s last two indices, producing the R-symbol that B encodes. Naively + # using a single `add_braid!` is wrong: it would resolve cyclic moves as crossings and + # pick up spurious R-symbol factors. + + A_in_layout = (oindA == codA && cindA == domA) + if A_in_layout + A′ = A + else + A′ = TO.tensoralloc_add( + scalartype(A), A, (oindA, cindA), false, Val(true), allocator + ) + add_transpose!(A′, A, (oindA, cindA), One(), Zero(), backend) end - p = (oindA, reverse(cindA)) - N = length(oindA) - levels = (ntuple(identity, N)..., (B.adjoint ? (N + 1, N + 2) : (N + 2, N + 1))...) - return add_braid!(C, A, p, levels, α, β, backend) + levelsB = B.adjoint ? (1, 2, 2, 1) : (2, 1, 1, 2) + N = numind(A) + M = N - 2 + levels = ( + ntuple(Returns(3), M)..., + levelsB[cindB[1]], levelsB[cindB[2]], + ) + + add_braid!( + C, A′, (ntuple(identity, M), (N, N - 1)), + levels, α, β, backend, + ) + + A_in_layout || TO.tensorfree!(A′, allocator) + return C end # ambiguity fix: @@ -325,277 +317,6 @@ function planarcontract!( α::Number, β::Number, backend, allocator ) return planarcontract!( - C, TensorMap(A), pA, TensorMap(B), pB, pAB, α, β, backend, allocator + C, A, pA, TensorMap(B), pB, pAB, α, β, backend, allocator ) end - -function planartrace!( - C::AbstractTensorMap, - A::BraidingTensor, - p::Index2Tuple, q::Index2Tuple, - α::Number, β::Number, - backend, allocator - ) - return planartrace!(C, TensorMap(A), p, q, α, β, backend, allocator) -end - -# function planarcontract!(C::AbstractTensorMap{<:Any,S,N₁,N₂}, -# A::BraidingTensor{S}, -# (oindA, cindA)::Index2Tuple{0,4}, -# B::AbstractTensorMap{<:Any,S}, -# (cindB, oindB)::Index2Tuple{4,<:Any}, -# (p1, p2)::Index2Tuple{N₁,N₂}, -# α::Number, β::Number, -# backend::Backend...) where {S,N₁,N₂} -# codA, domA = codomainind(A), domainind(A) -# codB, domB = codomainind(B), domainind(B) -# oindA, cindA, oindB, cindB = reorder_indices(codA, domA, codB, domB, oindA, cindA, -# oindB, cindB, p1, p2) - -# @assert space(B, cindB[1]) == space(A, cindA[1])' && -# space(B, cindB[2]) == space(A, cindA[2])' && -# space(B, cindB[3]) == space(A, cindA[3])' && -# space(B, cindB[4]) == space(A, cindA[4])' - -# if BraidingStyle(sectortype(B)) isa Bosonic -# return trace!(α, B, β, C, (), oindB, (cindB[1], cindB[2]), (cindB[3], cindB[4])) -# end - -# if iszero(β) -# fill!(C, β) -# elseif β != 1 -# rmul!(C, β) -# end -# I = sectortype(B) -# u = unit(I) -# f₀ = FusionTree{I}((), u, (), (), ()) -# braidingtensor_levels = A.adjoint ? (1, 2, 2, 1) : (2, 1, 1, 2) -# inv_braid = braidingtensor_levels[cindA[2]] > braidingtensor_levels[cindA[3]] -# for (f₁, f₂) in fusiontrees(B) -# local newtrees -# for ((f₁′, f₂′), coeff′) in transpose(f₁, f₂, cindB, oindB) -# f₁′.coupled == u || continue -# a = f₁′.uncoupled[1] -# b = f₁′.uncoupled[2] -# f₁′.uncoupled[3] == dual(a) || continue -# f₁′.uncoupled[4] == dual(b) || continue -# # should be automatic by matching spaces: -# # f₁′.isdual[1] != f₁′.isdual[3] || continue -# # f₁′.isdual[2] != f₁′.isdual[4] || continue -# for (f₁′′, coeff′′) in artin_braid(f₁′, 2; inv=inv_braid) -# f₁′′.innerlines[1] == u || continue -# coeff = coeff′ * coeff′′ * sqrtdim(a) * sqrtdim(b) -# if f₁′′.isdual[1] -# coeff *= frobenius_schur_phase(a) -# end -# if f₁′′.isdual[3] -# coeff *= frobenius_schur_phase(b) -# end -# f12 = (f₀, f₂′) -# if @isdefined newtrees -# newtrees[f12] = get(newtrees, f12, zero(coeff)) + coeff -# else -# newtrees = Dict(f12 => coeff) -# end -# end -# end -# @isdefined(newtrees) || continue -# for ((f₁′, f₂′), coeff) in newtrees -# TO._trace!(coeff * α, B[f₁, f₂], true, C[f₁′, f₂′], oindB, -# (cindB[1], cindB[2]), (cindB[3], cindB[4])) -# end -# end -# return C -# end -# function planarcontract!(C::AbstractTensorMap{<:Any,S,N₁,N₂}, -# A::AbstractTensorMap{<:Any,S}, -# (oindA, cindA)::Index2Tuple{0,4}, -# B::BraidingTensor{S}, -# (cindB, oindB)::Index2Tuple{4,<:Any}, -# (p1, p2)::Index2Tuple{N₁,N₂}, -# α::Number, β::Number, -# backends...) where {S,N₁,N₂} -# codA, domA = codomainind(A), domainind(A) -# codB, domB = codomainind(B), domainind(B) -# oindA, cindA, oindB, cindB = reorder_indices(codA, domA, codB, domB, oindA, cindA, -# oindB, cindB, p1, p2) - -# @assert space(B, cindB[1]) == space(A, cindA[1])' && -# space(B, cindB[2]) == space(A, cindA[2])' && -# space(B, cindB[3]) == space(A, cindA[3])' && -# space(B, cindB[4]) == space(A, cindA[4])' - -# if BraidingStyle(sectortype(B)) isa Bosonic -# return trace!(α, A, β, C, oindA, (), (cindA[1], cindA[2]), (cindA[3], cindA[4])) -# end - -# if iszero(β) -# fill!(C, β) -# elseif β != 1 -# rmul!(C, β) -# end -# I = sectortype(B) -# u = unit(I) -# f₀ = FusionTree{I}((), u, (), (), ()) -# braidingtensor_levels = B.adjoint ? (1, 2, 2, 1) : (2, 1, 1, 2) -# inv_braid = braidingtensor_levels[cindB[2]] > braidingtensor_levels[cindB[3]] -# for (f₁, f₂) in fusiontrees(A) -# local newtrees -# for ((f₁′, f₂′), coeff′) in transpose(f₁, f₂, oindA, cindA) -# f₂′.coupled == u || continue -# a = f₂′.uncoupled[1] -# b = f₂′.uncoupled[2] -# f₂′.uncoupled[3] == dual(a) || continue -# f₂′.uncoupled[4] == dual(b) || continue -# # should be automatic by matching spaces: -# # f₂′.isdual[1] != f₂′.isdual[3] || continue -# # f₂′.isdual[3] != f₂′.isdual[4] || continue -# for (f₂′′, coeff′′) in artin_braid(f₂′, 2; inv=inv_braid) -# f₂′′.innerlines[1] == u || continue -# coeff = coeff′ * conj(coeff′′ * sqrtdim(a) * sqrtdim(b)) -# if f₂′′.isdual[1] -# coeff *= conj(frobenius_schur_phase(a)) -# end -# if f₂′′.isdual[3] -# coeff *= conj(frobenius_schur_phase(b)) -# end -# f12 = (f₁′, f₀) -# if @isdefined newtrees -# newtrees[f12] = get(newtrees, f12, zero(coeff)) + coeff -# else -# newtrees = Dict(f12 => coeff) -# end -# end -# end -# @isdefined(newtrees) || continue -# for ((f₁′, f₂′), coeff) in newtrees -# TO._trace!(coeff * α, A[f₁, f₂], true, C[f₁′, f₂′], oindA, -# (cindA[1], cindA[2]), (cindA[3], cindA[4])) -# end -# end -# return C -# end -# function planarcontract!(C::AbstractTensorMap{<:Any,S,N₁,N₂}, -# A::BraidingTensor{S}, -# (oindA, cindA)::Index2Tuple{1,3}, -# B::AbstractTensorMap{<:Any,S}, -# (cindB, oindB)::Index2Tuple{1,<:Any}, -# (p1, p2)::Index2Tuple{N₁,N₂}, -# α::Number, β::Number, -# backend::Backend...) where {S,N₁,N₂} -# codA, domA = codomainind(A), domainind(A) -# codB, domB = codomainind(B), domainind(B) -# oindA, cindA, oindB, cindB = reorder_indices(codA, domA, codB, domB, oindA, cindA, -# oindB, cindB, p1, p2) - -# @assert space(B, cindB[1]) == space(A, cindA[1])' && -# space(B, cindB[2]) == space(A, cindA[2])' && -# space(B, cindB[3]) == space(A, cindA[3])' - -# if BraidingStyle(sectortype(B)) isa Bosonic -# return trace!(α, B, β, C, (cindB[2],), oindB, (cindB[1],), (cindB[3],)) -# end - -# if iszero(β) -# fill!(C, β) -# elseif β != 1 -# rmul!(C, β) -# end -# I = sectortype(B) -# u = unit(I) -# braidingtensor_levels = A.adjoint ? (1, 2, 2, 1) : (2, 1, 1, 2) -# inv_braid = braidingtensor_levels[cindA[2]] > braidingtensor_levels[cindA[3]] -# for (f₁, f₂) in fusiontrees(B) -# local newtrees -# for ((f₁′, f₂′), coeff′) in transpose(f₁, f₂, cindB, oindB) -# a = f₁′.uncoupled[1] -# b = f₁′.uncoupled[2] -# b == f₁′.coupled || continue -# a == dual(f₁′.uncoupled[3]) || continue -# # should be automatic by matching spaces: -# # f₁′.isdual[1] != f₁.isdual[3] || continue -# for (f₁′′, coeff′′) in artin_braid(f₁′, 2; inv=inv_braid) -# f₁′′.innerlines[1] == u || continue -# coeff = coeff′ * coeff′′ * sqrtdim(a) -# if f₁′′.isdual[1] -# coeff *= frobenius_schur_phase(a) -# end -# f₁′′′ = FusionTree{I}((b,), b, (f₁′′.isdual[3],), (), ()) -# f12 = (f₁′′′, f₂′) -# if @isdefined newtrees -# newtrees[f12] = get(newtrees, f12, zero(coeff)) + coeff -# else -# newtrees = Dict(f12 => coeff) -# end -# end -# end -# @isdefined(newtrees) || continue -# for ((f₁′, f₂′), coeff) in newtrees -# TO._trace!(coeff * α, B[f₁, f₂], true, C[f₁′, f₂′], -# (cindB[2], oindB...), (cindB[1],), (cindB[3],)) -# end -# end -# return C -# end -# function planarcontract!(C::AbstractTensorMap{<:Any,S,N₁,N₂}, -# A::AbstractTensorMap{<:Any,S}, -# (oindA, cindA)::Index2Tuple{<:Any,3}, -# B::BraidingTensor{S}, -# (cindB, oindB)::Index2Tuple{3,1}, -# (p1, p2)::Index2Tuple{N₁,N₂}, -# α::Number, β::Number, -# backend::Backend...) where {S,N₁,N₂} -# codA, domA = codomainind(A), domainind(A) -# codB, domB = codomainind(B), domainind(B) -# oindA, cindA, oindB, cindB = reorder_indices(codA, domA, codB, domB, oindA, cindA, -# oindB, cindB, p1, p2) - -# @assert space(B, cindB[1]) == space(A, cindA[1])' && -# space(B, cindB[2]) == space(A, cindA[2])' && -# space(B, cindB[3]) == space(A, cindA[3])' - -# if BraidingStyle(sectortype(A)) isa Bosonic -# return trace!(α, A, β, C, oindA, (cindA[2],), (cindA[1],), (cindA[3],)) -# end - -# if iszero(β) -# fill!(C, β) -# elseif β != 1 -# rmul!(C, β) -# end -# I = sectortype(B) -# u = unit(I) -# braidingtensor_levels = B.adjoint ? (1, 2, 2, 1) : (2, 1, 1, 2) -# inv_braid = braidingtensor_levels[cindB[2]] > braidingtensor_levels[cindB[3]] -# for (f₁, f₂) in fusiontrees(A) -# local newtrees -# for ((f₁′, f₂′), coeff′) in transpose(f₁, f₂, oindA, cindA) -# a = f₂′.uncoupled[1] -# b = f₂′.uncoupled[2] -# b == f₂′.coupled || continue -# a == dual(f₂′.uncoupled[3]) || continue -# # should be automatic by matching spaces: -# # f₂′.isdual[1] != f₂.isdual[3] || continue -# for (f₂′′, coeff′′) in artin_braid(f₂′, 2; inv=inv_braid) -# f₂′′.innerlines[1] == u || continue -# coeff = coeff′ * conj(coeff′′ * sqrtdim(a)) -# if f₂′′.isdual[1] -# coeff *= conj(frobenius_schur_phase(a)) -# end -# f₂′′′ = FusionTree{I}((b,), b, (f₂′′.isdual[3],), (), ()) -# f12 = (f₁′, f₂′′′) -# if @isdefined newtrees -# newtrees[f12] = get(newtrees, f12, zero(coeff)) + coeff -# else -# newtrees = Dict(f12 => coeff) -# end -# end -# end -# @isdefined(newtrees) || continue -# for ((f₁′, f₂′), coeff) in newtrees -# TO._trace!(coeff * α, A[f₁, f₂], true, C[f₁′, f₂′], -# (oindA..., cindA[2]), (cindA[1],), (cindA[3],)) -# end -# end -# return C -# end diff --git a/test/braidingtensor.jl b/test/braidingtensor.jl deleted file mode 100644 index 9b06b3536..000000000 --- a/test/braidingtensor.jl +++ /dev/null @@ -1,135 +0,0 @@ -# TODO: Make into proper tests and integrate in testset -# note: this is not part of the testsuite! - -import TensorKit: BraidingTensor - -V1 = GradedSpace{FermionSpin}(0 => 2, 1 / 2 => 2, 1 => 1, 3 / 2 => 1) - -V2 = GradedSpace{FibonacciAnyon}(:I => 2, :τ => 2) - -V3 = GradedSpace{IsingAnyon}(:I => 2, :psi => 1, :sigma => 1) - -for V in (V1, V2, V3) - @show V - - t = randn(V * V' * V' * V, V * V') - - ττ = TensorMap(BraidingTensor(V, V')) - @planar2 t1[-1 -2 -3 -4; -5 -6] := τ[-1 -2; 1 2] * t[1 2 -3 -4; -5 -6] - @planar2 t2[-1 -2 -3 -4; -5 -6] := ττ[-1 -2; 1 2] * t[1 2 -3 -4; -5 -6] - @planar2 t3[-1 -2 -3 -4; -5 -6] := τ[2 1; -2 -1] * t[1 2 -3 -4; -5 -6] - @planar2 t4[-1 -2 -3 -4; -5 -6] := τ'[-2 2; -1 1] * t[1 2 -3 -4; -5 -6] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V', V')) - @planar2 t1[-1 -2 -3 -4; -5 -6] := τ[-2 -3; 1 2] * t[-1 1 2 -4; -5 -6] - @planar2 t2[-1 -2 -3 -4; -5 -6] := ττ[-2 -3; 1 2] * t[-1 1 2 -4; -5 -6] - @planar2 t3[-1 -2 -3 -4; -5 -6] := τ[2 1; -3 -2] * t[-1 1 2 -4; -5 -6] - @planar2 t4[-1 -2 -3 -4; -5 -6] := τ'[-3 2; -2 1] * t[-1 1 2 -4; -5 -6] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V', V)) - @planar2 t1[-1 -2 -3 -4; -5 -6] := τ[1 -2; 2 -3] * t[-1 1 2 -4; -5 -6] - @planar2 t2[-1 -2 -3 -4; -5 -6] := ττ[1 -2; 2 -3] * t[-1 1 2 -4; -5 -6] - @planar2 t3[-1 -2 -3 -4; -5 -6] := τ[-3 2; -2 1] * t[-1 1 2 -4; -5 -6] - @planar2 t4[-1 -2 -3 -4; -5 -6] := τ'[-2 -3; 1 2] * t[-1 1 2 -4; -5 -6] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V, V')) - @planar2 t1[-1 -2 -3 -4; -5 -6] := τ[-3 2; -2 1] * t[-1 1 2 -4; -5 -6] - @planar2 t2[-1 -2 -3 -4; -5 -6] := ττ[-3 2; -2 1] * t[-1 1 2 -4; -5 -6] - @planar2 t3[-1 -2 -3 -4; -5 -6] := τ[1 -2; 2 -3] * t[-1 1 2 -4; -5 -6] - @planar2 t4[-1 -2 -3 -4; -5 -6] := τ'[2 1; -3 -2] * t[-1 1 2 -4; -5 -6] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V, V)) - @planar2 t1[-1 -2 -3 -4; -5 -6] := τ[2 1; -3 -2] * t[-1 1 2 -4; -5 -6] - @planar2 t2[-1 -2 -3 -4; -5 -6] := ττ[2 1; -3 -2] * t[-1 1 2 -4; -5 -6] - @planar2 t3[-1 -2 -3 -4; -5 -6] := τ[-2 -3; 1 2] * t[-1 1 2 -4; -5 -6] - @planar2 t4[-1 -2 -3 -4; -5 -6] := τ'[1 -2; 2 -3] * t[-1 1 2 -4; -5 -6] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V', V)) - @planar2 t1[-1 -2 -3 -4; -5 -6] := τ[-3 -4; 1 2] * t[-1 -2 1 2; -5 -6] - @planar2 t2[-1 -2 -3 -4; -5 -6] := ττ[-3 -4; 1 2] * t[-1 -2 1 2; -5 -6] - @planar2 t3[-1 -2 -3 -4; -5 -6] := τ[2 1; -4 -3] * t[-1 -2 1 2; -5 -6] - @planar2 t4[-1 -2 -3 -4; -5 -6] := τ'[-4 2; -3 1] * t[-1 -2 1 2; -5 -6] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V', V)) - @planar2 t1[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ[1 2; -5 -6] - @planar2 t2[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * ττ[1 2; -5 -6] - @planar2 t3[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ[-6 -5; 2 1] - @planar2 t4[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ'[2 -6; 1 -5] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V, V')) - @planar2 t1[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ'[1 2; -5 -6] - @planar2 t2[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * ττ'[1 2; -5 -6] - @planar2 t3[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ'[-6 -5; 2 1] - @planar2 t4[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ[2 -6; 1 -5] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V, V)) - @planar2 t1[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 1; -5 2] * τ[-4 -6; 1 2] - @planar2 t2[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 1; -5 2] * ττ[-4 -6; 1 2] - @planar2 t3[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 1; -5 2] * τ[2 1; -6 -4] - @planar2 t4[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 1; -5 2] * τ'[-6 2; -4 1] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V', V)) - @planar2 t1[(); (-1, -2)] := τ[2 1; 3 4] * t[1 2 3 4; -1 -2] - @planar2 t2[(); (-1, -2)] := ττ[2 1; 3 4] * t[1 2 3 4; -1 -2] - @planar2 t3[(); (-1, -2)] := τ[4 3; 1 2] * t[1 2 3 4; -1 -2] - @planar2 t4[(); (-1, -2)] := τ'[1 4; 2 3] * t[1 2 3 4; -1 -2] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V, V)) - @planar2 t1[-1; -2] := τ[2 1; 3 4] * t[-1 1 2 3; -2 4] - @planar2 t2[-1; -2] := ττ[2 1; 3 4] * t[-1 1 2 3; -2 4] - @planar2 t3[-1; -2] := τ[4 3; 1 2] * t[-1 1 2 3; -2 4] - @planar2 t4[-1; -2] := τ'[1 4; 2 3] * t[-1 1 2 3; -2 4] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V, V')) - @planar2 t1[-1 -2] := τ[2 1; 3 4] * t[-1 -2 1 2; 4 3] - @planar2 t2[-1 -2] := ττ[2 1; 3 4] * t[-1 -2 1 2; 4 3] - @planar2 t3[-1 -2] := τ[4 3; 1 2] * t[-1 -2 1 2; 4 3] - @planar2 t4[-1 -2] := τ'[1 4; 2 3] * t[-1 -2 1 2; 4 3] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V, V')) - @planar2 t1[-1 -2; -3 -4] := τ[-1 3; 1 2] * t[1 2 3 -2; -3 -4] - @planar2 t2[-1 -2; -3 -4] := ττ[-1 3; 1 2] * t[1 2 3 -2; -3 -4] - @planar2 t3[-1 -2; -3 -4] := τ[2 1; 3 -1] * t[1 2 3 -2; -3 -4] - @planar2 t4[-1 -2; -3 -4] := τ'[3 2; -1 1] * t[1 2 3 -2; -3 -4] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V', V')) - @planar2 t1[-1 -2; -3 -4] := τ'[-2 3; 1 2] * t[-1 1 2 3; -3 -4] - @planar2 t2[-1 -2; -3 -4] := ττ'[-2 3; 1 2] * t[-1 1 2 3; -3 -4] - @planar2 t3[-1 -2; -3 -4] := τ'[2 1; 3 -2] * t[-1 1 2 3; -3 -4] - @planar2 t4[-1 -2; -3 -4] := τ[3 2; -2 1] * t[-1 1 2 3; -3 -4] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V', V)) - @planar2 t1[-1 -2 -3; -4] := τ[-3 3; 1 2] * t[-1 -2 1 2; -4 3] - @planar2 t2[-1 -2 -3; -4] := ττ[-3 3; 1 2] * t[-1 -2 1 2; -4 3] - @planar2 t3[-1 -2 -3; -4] := τ[2 1; 3 -3] * t[-1 -2 1 2; -4 3] - @planar2 t4[-1 -2 -3; -4] := τ'[3 2; -3 1] * t[-1 -2 1 2; -4 3] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V', V)) - @planar2 t1[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ[1 2; -4 3] - @planar2 t2[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * ττ[1 2; -4 3] - @planar2 t3[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ[3 -4; 2 1] - @planar2 t4[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ'[2 3; 1 -4] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) - - ττ = TensorMap(BraidingTensor(V, V')) - @planar2 t1[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ'[1 2; -4 3] - @planar2 t2[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * ττ'[1 2; -4 3] - @planar2 t3[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ'[3 -4; 2 1] - @planar2 t4[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ[2 3; 1 -4] - @show norm(t1 - t2), norm(t1 - t3), norm(t1 - t4) -end diff --git a/test/runtests.jl b/test/runtests.jl index a7d2e2373..8bced2d1e 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -5,7 +5,6 @@ testsuite = ParallelTestRunner.find_tests(@__DIR__) # Exclude non-test files delete!(testsuite, "setup") # shared setup module -delete!(testsuite, "braidingtensor") # not part of the testsuite (see file header) # CUDA tests: only run if CUDA is functional using CUDA: CUDA diff --git a/test/tensors/braidingtensor.jl b/test/tensors/braidingtensor.jl new file mode 100644 index 000000000..bcf887b80 --- /dev/null +++ b/test/tensors/braidingtensor.jl @@ -0,0 +1,237 @@ +using Test, TestExtras +using TensorKit +using TensorKit: type_repr + +spacelist = default_spacelist(fast_tests) + +for V_tuple in spacelist + I = sectortype(first(V_tuple)) + Istr = type_repr(I) + BraidingStyle(I) isa NoBraiding && continue + V = first(V_tuple) + T = ComplexF64 + t = randn(T, V ⊗ V' ⊗ V' ⊗ V ← V ⊗ V') + println("---------------------------------------") + println("BraidingTensor planar contractions with symmetry: $Istr") + println("---------------------------------------") + @timedtestset "BraidingTensor planar contractions with symmetry: $Istr" verbose = true begin + @timedtestset "planaradd! with BraidingTensor" begin + b = BraidingTensor(V, V') + bb = TensorMap(b) + # Cyclic rotations of the planar leg cycle (cod1, cod2, dom2, dom1). + # Use transpose (F-symbols only) as reference, since permute requires SymmetricBraiding. + # rotation 0 (identity) + @planar t1[-1 -2; -3 -4] := b[-1 -2; -3 -4] + @test t1 ≈ bb + # rotation 1: single-tree cycle (4,1,2,3) → (p1=(2,4), p2=(1,3)) + @planar t2[-1 -2; -3 -4] := b[-3 -1; -4 -2] + @test t2 ≈ transpose(bb, ((2, 4), (1, 3))) + # rotation 2: single-tree cycle (3,4,1,2) → (p1=(4,3), p2=(2,1)) + @planar t3[-1 -2; -3 -4] := b[-4 -3; -2 -1] + @test t3 ≈ transpose(bb, ((4, 3), (2, 1))) + # rotation 3: single-tree cycle (2,3,4,1) → (p1=(3,1), p2=(4,2)) + @planar t4[-1 -2; -3 -4] := b[-2 -4; -1 -3] + @test t4 ≈ transpose(bb, ((3, 1), (4, 2))) + # adjoint BraidingTensor (rotation 0) + ba = b' + @planar t5[-1 -2; -3 -4] := ba[-1 -2; -3 -4] + @test t5 ≈ TensorMap(ba) + end + + @timedtestset "τ as left factor, all legs contracted" begin + # BraidingTensor(V, V') on leading codomain indices + ττ = TensorMap(BraidingTensor(V, V')) + @planar t1[-1 -2 -3 -4; -5 -6] := τ[-1 -2; 1 2] * t[1 2 -3 -4; -5 -6] + @planar t2[-1 -2 -3 -4; -5 -6] := ττ[-1 -2; 1 2] * t[1 2 -3 -4; -5 -6] + @planar t3[-1 -2 -3 -4; -5 -6] := τ[2 1; -2 -1] * t[1 2 -3 -4; -5 -6] + @planar t4[-1 -2 -3 -4; -5 -6] := τ'[-2 2; -1 1] * t[1 2 -3 -4; -5 -6] + @test t1 ≈ braid(t, ((2, 1, 3, 4), (5, 6)), (1, 2, 3, 4, 5, 6)) + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # BraidingTensor(V', V') on inner codomain indices + ττ = TensorMap(BraidingTensor(V', V')) + @planar t1[-1 -2 -3 -4; -5 -6] := τ[-2 -3; 1 2] * t[-1 1 2 -4; -5 -6] + @planar t2[-1 -2 -3 -4; -5 -6] := ττ[-2 -3; 1 2] * t[-1 1 2 -4; -5 -6] + @planar t3[-1 -2 -3 -4; -5 -6] := τ[2 1; -3 -2] * t[-1 1 2 -4; -5 -6] + @planar t4[-1 -2 -3 -4; -5 -6] := τ'[-3 2; -2 1] * t[-1 1 2 -4; -5 -6] + @test t1 ≈ braid(t, ((1, 3, 2, 4), (5, 6)), (1, 2, 3, 4, 5, 6)) + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # BraidingTensor(V', V) on trailing codomain indices + ττ = TensorMap(BraidingTensor(V', V)) + @planar t1[-1 -2 -3 -4; -5 -6] := τ[-3 -4; 1 2] * t[-1 -2 1 2; -5 -6] + @planar t2[-1 -2 -3 -4; -5 -6] := ττ[-3 -4; 1 2] * t[-1 -2 1 2; -5 -6] + @planar t3[-1 -2 -3 -4; -5 -6] := τ[2 1; -4 -3] * t[-1 -2 1 2; -5 -6] + @planar t4[-1 -2 -3 -4; -5 -6] := τ'[-4 2; -3 1] * t[-1 -2 1 2; -5 -6] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + end + + @timedtestset "τ as left factor, mixed open legs" begin + # BraidingTensor(V', V) with mixed index pattern + ττ = TensorMap(BraidingTensor(V', V)) + @planar t1[-1 -2 -3 -4; -5 -6] := τ[1 -2; 2 -3] * t[-1 1 2 -4; -5 -6] + @planar t2[-1 -2 -3 -4; -5 -6] := ττ[1 -2; 2 -3] * t[-1 1 2 -4; -5 -6] + @planar t3[-1 -2 -3 -4; -5 -6] := τ[-3 2; -2 1] * t[-1 1 2 -4; -5 -6] + @planar t4[-1 -2 -3 -4; -5 -6] := τ'[-2 -3; 1 2] * t[-1 1 2 -4; -5 -6] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # BraidingTensor(V, V') with mixed index pattern (inverse of previous) + ττ = TensorMap(BraidingTensor(V, V')) + @planar t1[-1 -2 -3 -4; -5 -6] := τ[-3 2; -2 1] * t[-1 1 2 -4; -5 -6] + @planar t2[-1 -2 -3 -4; -5 -6] := ττ[-3 2; -2 1] * t[-1 1 2 -4; -5 -6] + @planar t3[-1 -2 -3 -4; -5 -6] := τ[1 -2; 2 -3] * t[-1 1 2 -4; -5 -6] + @planar t4[-1 -2 -3 -4; -5 -6] := τ'[2 1; -3 -2] * t[-1 1 2 -4; -5 -6] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # BraidingTensor(V, V) with mixed index pattern + ττ = TensorMap(BraidingTensor(V, V)) + @planar t1[-1 -2 -3 -4; -5 -6] := τ[2 1; -3 -2] * t[-1 1 2 -4; -5 -6] + @planar t2[-1 -2 -3 -4; -5 -6] := ττ[2 1; -3 -2] * t[-1 1 2 -4; -5 -6] + @planar t3[-1 -2 -3 -4; -5 -6] := τ[-2 -3; 1 2] * t[-1 1 2 -4; -5 -6] + @planar t4[-1 -2 -3 -4; -5 -6] := τ'[1 -2; 2 -3] * t[-1 1 2 -4; -5 -6] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + end + + @timedtestset "τ as right factor" begin + # BraidingTensor(V', V) on all domain indices + ττ = TensorMap(BraidingTensor(V', V)) + @planar t1[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ[1 2; -5 -6] + @planar t2[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * ττ[1 2; -5 -6] + @planar t3[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ[-6 -5; 2 1] + @planar t4[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ'[2 -6; 1 -5] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # BraidingTensor(V, V') adjoint on all domain indices + ττ = TensorMap(BraidingTensor(V, V')) + @planar t1[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ'[1 2; -5 -6] + @planar t2[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * ττ'[1 2; -5 -6] + @planar t3[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ'[-6 -5; 2 1] + @planar t4[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 -4; 1 2] * τ[2 -6; 1 -5] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # BraidingTensor(V, V) with mixed domain legs + ττ = TensorMap(BraidingTensor(V, V)) + @planar t1[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 1; -5 2] * τ[-4 -6; 1 2] + @planar t2[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 1; -5 2] * ττ[-4 -6; 1 2] + @planar t3[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 1; -5 2] * τ[2 1; -6 -4] + @planar t4[-1 -2 -3 -4; -5 -6] := t[-1 -2 -3 1; -5 2] * τ'[-6 2; -4 1] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + end + + @timedtestset "τ with fully contracted output" begin + # scalar output + ττ = TensorMap(BraidingTensor(V', V)) + @planar t1[(); (-1, -2)] := τ[2 1; 3 4] * t[1 2 3 4; -1 -2] + @planar t2[(); (-1, -2)] := ττ[2 1; 3 4] * t[1 2 3 4; -1 -2] + @planar t3[(); (-1, -2)] := τ[4 3; 1 2] * t[1 2 3 4; -1 -2] + @planar t4[(); (-1, -2)] := τ'[1 4; 2 3] * t[1 2 3 4; -1 -2] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # rank-1 output + ττ = TensorMap(BraidingTensor(V, V)) + @planar t1[-1; -2] := τ[2 1; 3 4] * t[-1 1 2 3; -2 4] + @planar t2[-1; -2] := ττ[2 1; 3 4] * t[-1 1 2 3; -2 4] + @planar t3[-1; -2] := τ[4 3; 1 2] * t[-1 1 2 3; -2 4] + @planar t4[-1; -2] := τ'[1 4; 2 3] * t[-1 1 2 3; -2 4] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # rank-2 output + ττ = TensorMap(BraidingTensor(V, V')) + @planar t1[-1 -2] := τ[2 1; 3 4] * t[-1 -2 1 2; 4 3] + @planar t2[-1 -2] := ττ[2 1; 3 4] * t[-1 -2 1 2; 4 3] + @planar t3[-1 -2] := τ[4 3; 1 2] * t[-1 -2 1 2; 4 3] + @planar t4[-1 -2] := τ'[1 4; 2 3] * t[-1 -2 1 2; 4 3] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + end + + @timedtestset "τ with one open codomain leg" begin + # BraidingTensor(V, V') with one open codomain leg + ττ = TensorMap(BraidingTensor(V, V')) + @planar t1[-1 -2; -3 -4] := τ[-1 3; 1 2] * t[1 2 3 -2; -3 -4] + @planar t2[-1 -2; -3 -4] := ττ[-1 3; 1 2] * t[1 2 3 -2; -3 -4] + @planar t3[-1 -2; -3 -4] := τ[2 1; 3 -1] * t[1 2 3 -2; -3 -4] + @planar t4[-1 -2; -3 -4] := τ'[3 2; -1 1] * t[1 2 3 -2; -3 -4] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # BraidingTensor(V', V') adjoint with one open codomain leg + ττ = TensorMap(BraidingTensor(V', V')) + @planar t1[-1 -2; -3 -4] := τ'[-2 3; 1 2] * t[-1 1 2 3; -3 -4] + @planar t2[-1 -2; -3 -4] := ττ'[-2 3; 1 2] * t[-1 1 2 3; -3 -4] + @planar t3[-1 -2; -3 -4] := τ'[2 1; 3 -2] * t[-1 1 2 3; -3 -4] + @planar t4[-1 -2; -3 -4] := τ[3 2; -2 1] * t[-1 1 2 3; -3 -4] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # BraidingTensor(V', V) with one open codomain leg + ττ = TensorMap(BraidingTensor(V', V)) + @planar t1[-1 -2 -3; -4] := τ[-3 3; 1 2] * t[-1 -2 1 2; -4 3] + @planar t2[-1 -2 -3; -4] := ττ[-3 3; 1 2] * t[-1 -2 1 2; -4 3] + @planar t3[-1 -2 -3; -4] := τ[2 1; 3 -3] * t[-1 -2 1 2; -4 3] + @planar t4[-1 -2 -3; -4] := τ'[3 2; -3 1] * t[-1 -2 1 2; -4 3] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + end + + @timedtestset "τ as right factor with open domain leg" begin + # BraidingTensor(V', V) as right factor with one open domain leg + ττ = TensorMap(BraidingTensor(V', V)) + @planar t1[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ[1 2; -4 3] + @planar t2[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * ττ[1 2; -4 3] + @planar t3[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ[3 -4; 2 1] + @planar t4[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ'[2 3; 1 -4] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + + # BraidingTensor(V, V') adjoint as right factor with one open domain leg + ττ = TensorMap(BraidingTensor(V, V')) + @planar t1[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ'[1 2; -4 3] + @planar t2[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * ττ'[1 2; -4 3] + @planar t3[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ'[3 -4; 2 1] + @planar t4[-1 -2 -3; -4] := t[-1 -2 -3 3; 1 2] * τ[2 3; 1 -4] + @test t1 ≈ t2 + @test t1 ≈ t3 + @test t1 ≈ t4 + end + + @timedtestset "BraidingTensor × BraidingTensor" begin + # b1 domain == b2 codomain == V⊗V', straight-through (planar) contraction + b1 = BraidingTensor(V, V') # space: V'⊗V ← V⊗V' + b2 = BraidingTensor(V', V) # space: V⊗V' ← V'⊗V + bb1 = TensorMap(b1) + bb2 = TensorMap(b2) + @planar t1[-1 -2; -3 -4] := b1[-1 -2; 1 2] * b2[1 2; -3 -4] + @planar t2[-1 -2; -3 -4] := bb1[-1 -2; 1 2] * bb2[1 2; -3 -4] + @test t1 ≈ t2 + end + end + TensorKit.empty_globalcaches!() +end