diff --git a/src/sectors/fermions.jl b/src/sectors/fermions.jl index e78ea498..bd6d78f7 100644 --- a/src/sectors/fermions.jl +++ b/src/sectors/fermions.jl @@ -1,57 +1,51 @@ -struct Fermion{P,I<:Sector} <: Sector - sector::I - function Fermion{P,I}(sector::I) where {P, I<:Sector} - @assert BraidingStyle(I) isa Bosonic - return new{P,I}(sector) - end +struct FermionParity <: Sector + sector::Z2Irrep end -Fermion{P}(sector::I) where {P, I<:Sector} = Fermion{P,I}(sector) -Fermion{P,I}(sector) where {P, I<:Sector} = Fermion{P,I}(convert(I, sector)) -Base.convert(::Type{Fermion{P,I}}, a::Fermion{P,I}) where {P, I<:Sector} = a -Base.convert(::Type{Fermion{P,I}}, a) where {P, I<:Sector} = Fermion{P,I}(convert(I, a)) +const fℤ₂ = FermionParity +fermionparity(f::FermionParity) = isodd(f.sector.n) -fermionparity(f::Fermion{P}) where P = P(f.sector) +Base.convert(::Type{FermionParity}, a::FermionParity) = a +Base.convert(::Type{FermionParity}, a) = FermionParity(a) -Base.IteratorSize(::Type{SectorValues{Fermion{P,I}}}) where {P, I<:Sector} = - Base.IteratorSize(SectorValues{I}) -Base.length(::SectorValues{Fermion{P,I}}) where {P, I<:Sector} = length(values(I)) -function Base.iterate(::SectorValues{Fermion{P, I}}) where {P, I<:Sector} - next = iterate(values(I)) +Base.IteratorSize(::Type{SectorValues{FermionParity}}) = + Base.IteratorSize(SectorValues{Z2Irrep}) +Base.length(::SectorValues{FermionParity}) = length(values(Z2Irrep)) +function Base.iterate(::SectorValues{FermionParity}) + next = iterate(values(Z2Irrep)) @assert next !== nothing value, state = next - return Fermion{P}(value), state + return FermionParity(value), state end -function Base.iterate(::SectorValues{Fermion{P, I}}, state) where {P, I<:Sector} - next = iterate(values(I), state) + +function Base.iterate(::SectorValues{FermionParity}, state) + next = iterate(values(Z2Irrep), state) if next === nothing return nothing else value, state = next - return Fermion{P}(value), state + return FermionParity(value), state end end -Base.getindex(::SectorValues{Fermion{P, I}}, i) where {P, I<:Sector} = - Fermion{P}(values(I)[i]) -findindex(::SectorValues{Fermion{P, I}}, f::Fermion{P, I}) where {P, I<:Sector} = - findindex(values(I), f.sector) +Base.getindex(::SectorValues{FermionParity}, i) = FermionParity(values(Z2Irrep)[i]) +findindex(::SectorValues{FermionParity}, f::FermionParity) = findindex(values(Z2Irrep), f.sector) -Base.one(::Type{Fermion{P, I}}) where {P, I<:Sector} = Fermion{P}(one(I)) -Base.conj(f::Fermion{P}) where {P} = Fermion{P}(conj(f.sector)) +Base.one(::Type{FermionParity}) = FermionParity(one(Z2Irrep)) +Base.conj(f::FermionParity) = FermionParity(conj(f.sector)) -dim(f::Fermion) = dim(f.sector) +dim(f::FermionParity) = dim(f.sector) -FusionStyle(::Type{<:Fermion{<:Any,I}}) where {I<:Sector} = FusionStyle(I) -BraidingStyle(::Type{<:Fermion}) = Fermionic() -Base.isreal(::Type{Fermion{<:Any,I}}) where {I<:Sector} = isreal(I) +FusionStyle(::Type{FermionParity}) = FusionStyle(Z2Irrep) +BraidingStyle(::Type{FermionParity}) = Fermionic() +Base.isreal(::Type{FermionParity}) = isreal(Z2Irrep) -⊗(a::F, b::F) where {F<:Fermion} = SectorSet{F}(a.sector ⊗ b.sector) +⊗(a::FermionParity, b::FermionParity) = SectorSet{FermionParity}(a.sector ⊗ b.sector) -Nsymbol(a::F, b::F, c::F) where {F<:Fermion} = Nsymbol(a.sector, b.sector, c.sector) +Nsymbol(a::FermionParity, b::FermionParity, c::FermionParity) = Nsymbol(a.sector, b.sector, c.sector) -Fsymbol(a::F, b::F, c::F, d::F, e::F, f::F) where {F<:Fermion} = +Fsymbol(a::F, b::F, c::F, d::F, e::F, f::F) where {F<:FermionParity} = Fsymbol(a.sector, b.sector, c.sector, d.sector, e.sector, f.sector) -function Rsymbol(a::F, b::F, c::F) where {F<:Fermion} +function Rsymbol(a::F, b::F, c::F) where {F<:FermionParity} if fermionparity(a) && fermionparity(b) return -Rsymbol(a.sector, b.sector, c.sector) else @@ -59,34 +53,19 @@ function Rsymbol(a::F, b::F, c::F) where {F<:Fermion} end end -twist(a::Fermion) = ifelse(fermionparity(a), -1, +1)*twist(a.sector) +twist(a::FermionParity) = ifelse(fermionparity(a), -1, +1) * twist(a.sector) -type_repr(::Type{Fermion{P,I}}) where {P, I<:Sector} = "Fermion{$P, " * type_repr(I) * "}" +type_repr(::Type{FermionParity}) = "FermionParity" -function Base.show(io::IO, a::Fermion{P, I}) where {P, I<:Sector} - if get(io, :typeinfo, nothing) !== Fermion{P, I} +function Base.show(io::IO, a::FermionParity) + if get(io, :typeinfo, nothing) !== FermionParity print(io, type_repr(typeof(a)), "(") end - print(IOContext(io, :typeinfo => I), a.sector) - if get(io, :typeinfo, nothing) !== Fermion{P, I} + print(IOContext(io, :typeinfo => Z2Irrep), a.sector) + if get(io, :typeinfo, nothing) !== FermionParity print(io, ")") end end -Base.hash(f::Fermion, h::UInt) = hash(f.sector, h) -Base.isless(a::F, b::F) where {F<:Fermion} = isless(a.sector, b.sector) - -_fermionparity(a::Z2Irrep) = isodd(a.n) -_fermionnumber(a::U1Irrep) = isodd(convert(Int, a.charge)) -_fermionspin(a::SU2Irrep) = isodd(twice(a.j)) - -const FermionParity = Fermion{_fermionparity, Z2Irrep} -const FermionNumber = Fermion{_fermionnumber, U1Irrep} -const FermionSpin = Fermion{_fermionspin, SU2Irrep} -const fℤ₂ = FermionParity -const fU₁ = FermionNumber -const fSU₂ = FermionSpin - -type_repr(::Type{FermionParity}) = "FermionParity" -type_repr(::Type{FermionNumber}) = "FermionNumber" -type_repr(::Type{FermionSpin}) = "FermionSpin" +Base.hash(f::FermionParity, h::UInt) = hash(f.sector, h) +Base.isless(a::FermionParity, b::FermionParity) = isless(a.sector, b.sector) diff --git a/test/runtests.jl b/test/runtests.jl index 380f8aa5..b2a2aa71 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -19,7 +19,6 @@ const TK = TensorKit Random.seed!(1234) smallset(::Type{I}) where {I<:Sector} = take(values(I), 5) -smallset(::Type{FermionNumber}) = FermionNumber.((0, +1, -1, +2, -2)) function smallset(::Type{ProductSector{Tuple{I1,I2}}}) where {I1,I2} iter = product(smallset(I1), smallset(I2)) s = collect(i ⊠ j for (i,j) in iter if dim(i)*dim(j) <= 6) @@ -52,10 +51,11 @@ function hasfusiontensor(I::Type{<:Sector}) end sectorlist = (Z2Irrep, Z3Irrep, Z4Irrep, U1Irrep, CU1Irrep, SU2Irrep, NewSU2Irrep, SU3Irrep, - FibonacciAnyon, IsingAnyon, FermionParity, FermionNumber, FermionSpin, - FermionParity ⊠ FermionParity, Z3Irrep ⊠ Z4Irrep, FermionNumber ⊠ SU2Irrep, - FermionSpin ⊠ SU2Irrep, NewSU2Irrep ⊠ NewSU2Irrep, NewSU2Irrep ⊠ SU2Irrep, - FermionSpin ⊠ NewSU2Irrep, Z2Irrep ⊠ FibonacciAnyon ⊠ FibonacciAnyon) + FibonacciAnyon, IsingAnyon, FermionParity, FermionParity ⊠ FermionParity, + Z3Irrep ⊠ Z4Irrep, FermionParity ⊠ U1Irrep ⊠ SU2Irrep, + FermionParity ⊠ SU2Irrep ⊠ SU2Irrep, NewSU2Irrep ⊠ NewSU2Irrep, + NewSU2Irrep ⊠ SU2Irrep, FermionParity ⊠ SU2Irrep ⊠ NewSU2Irrep, + Z2Irrep ⊠ FibonacciAnyon ⊠ FibonacciAnyon) Ti = time() include("sectors.jl") diff --git a/test/tensors.jl b/test/tensors.jl index 882813f7..2c103220 100644 --- a/test/tensors.jl +++ b/test/tensors.jl @@ -1,53 +1,53 @@ Vtr = (ℂ^3, - (ℂ^4)', - ℂ^5, - ℂ^6, - (ℂ^7)') -Vℤ₂ = (ℂ[Z2Irrep](0=>1, 1=>1), - ℂ[Z2Irrep](0=>1, 1=>2)', - ℂ[Z2Irrep](0=>3, 1=>2)', - ℂ[Z2Irrep](0=>2, 1=>3), - ℂ[Z2Irrep](0=>2, 1=>5)) -Vfℤ₂ = (ℂ[FermionParity](0=>1, 1=>1), - ℂ[FermionParity](0=>1, 1=>2)', - ℂ[FermionParity](0=>3, 1=>2)', - ℂ[FermionParity](0=>2, 1=>3), - ℂ[FermionParity](0=>2, 1=>5)) -Vℤ₃ = (ℂ[Z3Irrep](0=>1, 1=>2, 2=>2), - ℂ[Z3Irrep](0=>3, 1=>1, 2=>1), - ℂ[Z3Irrep](0=>2, 1=>2, 2=>1)', - ℂ[Z3Irrep](0=>1, 1=>2, 2=>3), - ℂ[Z3Irrep](0=>1, 1=>3, 2=>3)') -VU₁ = (ℂ[U1Irrep](0=>1, 1=>2, -1=>2), - ℂ[U1Irrep](0=>3, 1=>1, -1=>1), - ℂ[U1Irrep](0=>2, 1=>2, -1=>1)', - ℂ[U1Irrep](0=>1, 1=>2, -1=>3), - ℂ[U1Irrep](0=>1, 1=>3, -1=>3)') -VfU₁ = (ℂ[FermionNumber](0=>1, 1=>2, -1=>2), - ℂ[FermionNumber](0=>3, 1=>1, -1=>1), - ℂ[FermionNumber](0=>2, 1=>2, -1=>1)', - ℂ[FermionNumber](0=>1, 1=>2, -1=>3), - ℂ[FermionNumber](0=>1, 1=>3, -1=>3)') -VCU₁ = (ℂ[CU1Irrep]((0,0)=>1, (0,1)=>2, 1=>1), - ℂ[CU1Irrep]((0,0)=>3, (0,1)=>0, 1=>1), - ℂ[CU1Irrep]((0,0)=>1, (0,1)=>0, 1=>2)', - ℂ[CU1Irrep]((0,0)=>2, (0,1)=>2, 1=>1), - ℂ[CU1Irrep]((0,0)=>2, (0,1)=>1, 1=>2)') -VSU₂ = (ℂ[SU2Irrep](0=>3, 1//2=>1), - ℂ[SU2Irrep](0=>2, 1=>1), - ℂ[SU2Irrep](1//2=>1, 1=>1)', - ℂ[SU2Irrep](0=>2, 1//2=>2), - ℂ[SU2Irrep](0=>1, 1//2=>1, 3//2=>1)') -VfSU₂ = (ℂ[FermionSpin](0=>3, 1//2=>1), - ℂ[FermionSpin](0=>2, 1=>1), - ℂ[FermionSpin](1//2=>1, 1=>1)', - ℂ[FermionSpin](0=>2, 1//2=>2), - ℂ[FermionSpin](0=>1, 1//2=>1, 3//2=>1)') -VSU₃ = (ℂ[SU3Irrep]((0,0,0)=>3, (1,0,0)=>1), - ℂ[SU3Irrep]((0,0,0)=>3, (2,0,0)=>1)', - ℂ[SU3Irrep]((1,1,0)=>1, (2,1,0)=>1), - ℂ[SU3Irrep]((1,0,0)=>1, (2,0,0)=>1), - ℂ[SU3Irrep]((0,0,0)=>1, (1,0,0)=>1, (1,1,0)=>1)') + (ℂ^4)', + ℂ^5, + ℂ^6, + (ℂ^7)') +Vℤ₂ = (ℂ[Z2Irrep](0 => 1, 1 => 1), + ℂ[Z2Irrep](0 => 1, 1 => 2)', + ℂ[Z2Irrep](0 => 3, 1 => 2)', + ℂ[Z2Irrep](0 => 2, 1 => 3), + ℂ[Z2Irrep](0 => 2, 1 => 5)) +Vfℤ₂ = (ℂ[FermionParity](0 => 1, 1 => 1), + ℂ[FermionParity](0 => 1, 1 => 2)', + ℂ[FermionParity](0 => 3, 1 => 2)', + ℂ[FermionParity](0 => 2, 1 => 3), + ℂ[FermionParity](0 => 2, 1 => 5)) +Vℤ₃ = (ℂ[Z3Irrep](0 => 1, 1 => 2, 2 => 2), + ℂ[Z3Irrep](0 => 3, 1 => 1, 2 => 1), + ℂ[Z3Irrep](0 => 2, 1 => 2, 2 => 1)', + ℂ[Z3Irrep](0 => 1, 1 => 2, 2 => 3), + ℂ[Z3Irrep](0 => 1, 1 => 3, 2 => 3)') +VU₁ = (ℂ[U1Irrep](0 => 1, 1 => 2, -1 => 2), + ℂ[U1Irrep](0 => 3, 1 => 1, -1 => 1), + ℂ[U1Irrep](0 => 2, 1 => 2, -1 => 1)', + ℂ[U1Irrep](0 => 1, 1 => 2, -1 => 3), + ℂ[U1Irrep](0 => 1, 1 => 3, -1 => 3)') +VfU₁ = (ℂ[FermionParity⊠U1Irrep]((0, 0) => 1, (1, 1) => 2, (1, -1) => 2), + ℂ[FermionParity⊠U1Irrep]((0, 0) => 3, (1, 1) => 1, (1, -1) => 1), + ℂ[FermionParity⊠U1Irrep]((0, 0) => 2, (1, 1) => 2, (1, -1) => 1)', + ℂ[FermionParity⊠U1Irrep]((0, 0) => 1, (1, 1) => 2, (1, -1) => 3), + ℂ[FermionParity⊠U1Irrep]((0, 0) => 1, (1, 1) => 3, (1, -1) => 3)') +VCU₁ = (ℂ[CU1Irrep]((0, 0) => 1, (0, 1) => 2, 1 => 1), + ℂ[CU1Irrep]((0, 0) => 3, (0, 1) => 0, 1 => 1), + ℂ[CU1Irrep]((0, 0) => 1, (0, 1) => 0, 1 => 2)', + ℂ[CU1Irrep]((0, 0) => 2, (0, 1) => 2, 1 => 1), + ℂ[CU1Irrep]((0, 0) => 2, (0, 1) => 1, 1 => 2)') +VSU₂ = (ℂ[SU2Irrep](0 => 3, 1 // 2 => 1), + ℂ[SU2Irrep](0 => 2, 1 => 1), + ℂ[SU2Irrep](1 // 2 => 1, 1 => 1)', + ℂ[SU2Irrep](0 => 2, 1 // 2 => 2), + ℂ[SU2Irrep](0 => 1, 1 // 2 => 1, 3 // 2 => 1)') +VfSU₂ = (ℂ[FermionParity⊠SU2Irrep]((0, 0) => 3, (1, 1 // 2) => 1), + ℂ[FermionParity⊠SU2Irrep]((0, 0) => 2, (0, 1) => 1), + ℂ[FermionParity⊠SU2Irrep]((1, 1 // 2) => 1, (0, 1) => 1)', + ℂ[FermionParity⊠SU2Irrep]((0, 0) => 2, (1, 1 // 2) => 2), + ℂ[FermionParity⊠SU2Irrep]((0, 0) => 1, (1, 1 // 2) => 1, (1, 3 // 2) => 1)') +VSU₃ = (ℂ[SU3Irrep]((0, 0, 0) => 3, (1, 0, 0) => 1), + ℂ[SU3Irrep]((0, 0, 0) => 3, (2, 0, 0) => 1)', + ℂ[SU3Irrep]((1, 1, 0) => 1, (2, 1, 0) => 1), + ℂ[SU3Irrep]((1, 0, 0) => 1, (2, 0, 0) => 1), + ℂ[SU3Irrep]((0, 0, 0) => 1, (1, 0, 0) => 1, (1, 1, 0) => 1)') for V in (Vtr, Vℤ₂, Vfℤ₂, Vℤ₃, VU₁, VfU₁, VCU₁, VSU₂, VfSU₂, VSU₃) V1, V2, V3, V4, V5 = V @@ -106,7 +106,7 @@ for V in spacelist W = V1 ⊗ V2 ⊗ V3 ← V4 ⊗ V5 for T in (Int, Float32, ComplexF64) if T == Int - t = TensorMap(sz->rand(-20:20, sz), W) + t = TensorMap(sz -> rand(-20:20, sz), W) else t = TensorMap(randn, T, W) end @@ -126,22 +126,22 @@ for V in spacelist @test codomain(t) == codomain(W) @test domain(t) == domain(W) @test isa(@constinferred(norm(t)), real(T)) - @test norm(t)^2 ≈ dot(t,t) + @test norm(t)^2 ≈ dot(t, t) α = rand(T) - @test norm(α*t) ≈ abs(α)*norm(t) - @test norm(t+t, 2) ≈ 2*norm(t, 2) - @test norm(t+t, 1) ≈ 2*norm(t, 1) - @test norm(t+t, Inf) ≈ 2*norm(t, Inf) - p = 3*rand(Float64) - @test norm(t+t, p) ≈ 2*norm(t, p) + @test norm(α * t) ≈ abs(α) * norm(t) + @test norm(t + t, 2) ≈ 2 * norm(t, 2) + @test norm(t + t, 1) ≈ 2 * norm(t, 1) + @test norm(t + t, Inf) ≈ 2 * norm(t, Inf) + p = 3 * rand(Float64) + @test norm(t + t, p) ≈ 2 * norm(t, p) @test norm(t) ≈ norm(t') t2 = TensorMap(rand, T, W) β = rand(T) - @test @constinferred(dot(β*t2,α*t)) ≈ conj(β)*α*conj(dot(t,t2)) - @test dot(t2,t) ≈ conj(dot(t, t2)) - @test dot(t2,t) ≈ conj(dot(t2', t')) - @test dot(t2,t) ≈ dot(t', t2') + @test @constinferred(dot(β * t2, α * t)) ≈ conj(β) * α * conj(dot(t, t2)) + @test dot(t2, t) ≈ conj(dot(t, t2)) + @test dot(t2, t) ≈ conj(dot(t2', t')) + @test dot(t2, t) ≈ dot(t', t2') i1 = @constinferred(isomorphism(Matrix{T}, V1 ⊗ V2, V2 ⊗ V1)) i2 = @constinferred(isomorphism(Matrix{T}, V2 ⊗ V1, V1 ⊗ V2)) @@ -150,9 +150,9 @@ for V in spacelist w = @constinferred(isometry(Matrix{T}, V1 ⊗ (oneunit(V1) ⊕ oneunit(V1)), V1)) - @test dim(w) == 2*dim(V1←V1) - @test w'*w == id(Matrix{T}, V1) - @test w*w' == (w*w')^2 + @test dim(w) == 2 * dim(V1 ← V1) + @test w' * w == id(Matrix{T}, V1) + @test w * w' == (w * w')^2 end end if hasfusiontensor(I) @@ -161,11 +161,11 @@ for V in spacelist for T in (Float32, ComplexF64) t = TensorMap(rand, T, W) t2 = TensorMap(rand, T, W) - @test norm(t, 2) ≈ norm(convert(Array,t), 2) - @test dot(t2,t) ≈ dot(convert(Array,t2), convert(Array, t)) + @test norm(t, 2) ≈ norm(convert(Array, t), 2) + @test dot(t2, t) ≈ dot(convert(Array, t2), convert(Array, t)) α = rand(T) - @test convert(Array, α*t) ≈ α*convert(Array,t) - @test convert(Array, t+t) ≈ 2*convert(Array,t) + @test convert(Array, α * t) ≈ α * convert(Array, t) + @test convert(Array, t + t) ≈ 2 * convert(Array, t) end end @timedtestset "Real and imaginary parts" begin @@ -192,12 +192,12 @@ for V in spacelist t′ = Tensor(rand, ComplexF64, W) for k = 0:5 for p in permutations(1:5) - p1 = ntuple(n->p[n], k) - p2 = ntuple(n->p[k+n], 5-k) + p1 = ntuple(n -> p[n], k) + p2 = ntuple(n -> p[k+n], 5 - k) t2 = @constinferred permute(t, p1, p2) @test norm(t2) ≈ norm(t) - t2′= permute(t′, p1, p2) - @test dot(t2′,t2) ≈ dot(t′,t) ≈ dot(transpose(t2′), transpose(t2)) + t2′ = permute(t′, p1, p2) + @test dot(t2′, t2) ≈ dot(t′, t) ≈ dot(transpose(t2′), transpose(t2)) end end end @@ -207,19 +207,19 @@ for V in spacelist t = Tensor(rand, ComplexF64, W) for k = 0:5 for p in permutations(1:5) - p1 = ntuple(n->p[n], k) - p2 = ntuple(n->p[k+n], 5-k) + p1 = ntuple(n -> p[n], k) + p2 = ntuple(n -> p[k+n], 5 - k) t2 = permute(t, p1, p2) a2 = convert(Array, t2) - @test a2 ≈ permutedims(convert(Array, t), (p1...,p2...)) - @test convert(Array, transpose(t2)) ≈ permutedims(a2, (5,4,3,2,1)) + @test a2 ≈ permutedims(convert(Array, t), (p1..., p2...)) + @test convert(Array, transpose(t2)) ≈ permutedims(a2, (5, 4, 3, 2, 1)) end end end end @timedtestset "Full trace: test self-consistency" begin t = Tensor(rand, ComplexF64, V1 ⊗ V2' ⊗ V2 ⊗ V1') - t2 = permute(t, (1,2), (4,3)) + t2 = permute(t, (1, 2), (4, 3)) s = @constinferred tr(t2) @test conj(s) ≈ tr(t2') if !isdual(V1) @@ -229,24 +229,24 @@ for V in spacelist t2 = twist!(t2, 2) end ss = tr(t2) - @tensor s2 = t[a,b,b,a] - @tensor t3[a,b] := t[a,c,c,b] - @tensor s3 = t3[a,a] + @tensor s2 = t[a, b, b, a] + @tensor t3[a, b] := t[a, c, c, b] + @tensor s3 = t3[a, a] @test ss ≈ s2 @test ss ≈ s3 end @timedtestset "Partial trace: test self-consistency" begin t = Tensor(rand, ComplexF64, V1 ⊗ V2' ⊗ V3 ⊗ V2 ⊗ V1' ⊗ V3') - @tensor t2[a,b] := t[c,d,b,d,c,a] - @tensor t4[a,b,c,d] := t[d,e,b,e,c,a] - @tensor t5[a,b] := t4[a,b,c,c] + @tensor t2[a, b] := t[c, d, b, d, c, a] + @tensor t4[a, b, c, d] := t[d, e, b, e, c, a] + @tensor t5[a, b] := t4[a, b, c, c] @test t2 ≈ t5 end if hasfusiontensor(I) @timedtestset "Trace: test via conversion" begin t = Tensor(rand, ComplexF64, V1 ⊗ V2' ⊗ V3 ⊗ V2 ⊗ V1' ⊗ V3') - @tensor t2[a,b] := t[c,d,b,d,c,a] - @tensor t3[a,b] := convert(Array, t)[c,d,b,d,c,a] + @tensor t2[a, b] := t[c, d, b, d, c, a] + @tensor t3[a, b] := convert(Array, t)[c, d, b, d, c, a] @test t3 ≈ convert(Array, t2) end end @@ -254,25 +254,25 @@ for V in spacelist t1 = Tensor(rand, ComplexF64, V1 ⊗ V2 ⊗ V3) t2 = Tensor(rand, ComplexF64, V2' ⊗ V4 ⊗ V1') t3 = t1 ⊗ t2 - @tensor ta[a,b] := t1[x,y,a]*t2[y,b,x] - @tensor tb[a,b] := t3[x,y,a,y,b,x] + @tensor ta[a, b] := t1[x, y, a] * t2[y, b, x] + @tensor tb[a, b] := t3[x, y, a, y, b, x] @test ta ≈ tb end if hasfusiontensor(I) @timedtestset "Tensor contraction: test via conversion" begin - A1 = TensorMap(randn, ComplexF64, V1'*V2', V3') - A2 = TensorMap(randn, ComplexF64, V3*V4, V5) + A1 = TensorMap(randn, ComplexF64, V1' * V2', V3') + A2 = TensorMap(randn, ComplexF64, V3 * V4, V5) rhoL = TensorMap(randn, ComplexF64, V1, V1) rhoR = TensorMap(randn, ComplexF64, V5, V5)' # test adjoint tensor - H = TensorMap(randn, ComplexF64, V2*V4, V2*V4) + H = TensorMap(randn, ComplexF64, V2 * V4, V2 * V4) @tensor HrA12[a, s1, s2, c] := rhoL[a, a'] * conj(A1[a', t1, b]) * - A2[b, t2, c'] * rhoR[c', c] * H[s1, s2, t1, t2] + A2[b, t2, c'] * rhoR[c', c] * H[s1, s2, t1, t2] @tensor HrA12array[a, s1, s2, c] := convert(Array, rhoL)[a, a'] * - conj(convert(Array, A1)[a', t1, b]) * - convert(Array, A2)[b, t2, c'] * - convert(Array, rhoR)[c', c] * - convert(Array, H)[s1, s2, t1, t2] + conj(convert(Array, A1)[a', t1, b]) * + convert(Array, A2)[b, t2, c'] * + convert(Array, rhoR)[c', c] * + convert(Array, H)[s1, s2, t1, t2] @test HrA12array ≈ convert(Array, HrA12) end @@ -284,15 +284,15 @@ for V in spacelist t1 = TensorMap(rand, T, W1, W1) t2 = TensorMap(rand, T, W2, W2) t = TensorMap(rand, T, W1, W2) - @test t1*(t1\t) ≈ t - @test (t/t2)*t2 ≈ t - @test t1\one(t1) ≈ inv(t1) - @test one(t1)/t1 ≈ pinv(t1) + @test t1 * (t1 \ t) ≈ t + @test (t / t2) * t2 ≈ t + @test t1 \ one(t1) ≈ inv(t1) + @test one(t1) / t1 ≈ pinv(t1) @test_throws SpaceMismatch inv(t) - @test_throws SpaceMismatch t2\t - @test_throws SpaceMismatch t/t1 - tp = pinv(t)*t - @test tp ≈ tp*tp + @test_throws SpaceMismatch t2 \ t + @test_throws SpaceMismatch t / t1 + tp = pinv(t) * t + @test tp ≈ tp * tp end end if hasfusiontensor(I) @@ -308,10 +308,10 @@ for V in spacelist At1 = reshape(convert(Array, t1), d1, d1) At2 = reshape(convert(Array, t2), d2, d2) At = reshape(convert(Array, t), d1, d2) - @test reshape(convert(Array, t1*t), d1, d2) ≈ At1*At - @test reshape(convert(Array, t1'*t), d1, d2) ≈ At1'*At - @test reshape(convert(Array, t2*t'), d2, d1) ≈ At2*At' - @test reshape(convert(Array, t2'*t'), d2, d1) ≈ At2'*At' + @test reshape(convert(Array, t1 * t), d1, d2) ≈ At1 * At + @test reshape(convert(Array, t1' * t), d1, d2) ≈ At1' * At + @test reshape(convert(Array, t2 * t'), d2, d1) ≈ At2 * At' + @test reshape(convert(Array, t2' * t'), d2, d1) ≈ At2' * At' @test reshape(convert(Array, inv(t1)), d1, d1) ≈ inv(At1) @test reshape(convert(Array, pinv(t)), d2, d1) ≈ pinv(At) @@ -320,15 +320,15 @@ for V in spacelist continue end - @test reshape(convert(Array, t1\t), d1, d2) ≈ At1\At - @test reshape(convert(Array, t1'\t), d1, d2) ≈ At1'\At - @test reshape(convert(Array, t2\t'), d2, d1) ≈ At2\At' - @test reshape(convert(Array, t2'\t'), d2, d1) ≈ At2'\At' + @test reshape(convert(Array, t1 \ t), d1, d2) ≈ At1 \ At + @test reshape(convert(Array, t1' \ t), d1, d2) ≈ At1' \ At + @test reshape(convert(Array, t2 \ t'), d2, d1) ≈ At2 \ At' + @test reshape(convert(Array, t2' \ t'), d2, d1) ≈ At2' \ At' - @test reshape(convert(Array, t2/t), d2, d1) ≈ At2/At - @test reshape(convert(Array, t2'/t), d2, d1) ≈ At2'/At - @test reshape(convert(Array, t1/t'), d1, d2) ≈ At1/At' - @test reshape(convert(Array, t1'/t'), d1, d2) ≈ At1'/At' + @test reshape(convert(Array, t2 / t), d2, d1) ≈ At2 / At + @test reshape(convert(Array, t2' / t), d2, d1) ≈ At2' / At + @test reshape(convert(Array, t1 / t'), d1, d2) ≈ At1 / At' + @test reshape(convert(Array, t1' / t'), d1, d2) ≈ At1' / At' end end end @@ -339,70 +339,70 @@ for V in spacelist ts = (Tensor(rand, T, W), Tensor(rand, T, W)') for t in ts @testset "leftorth with $alg" for alg in (TensorKit.QR(), TensorKit.QRpos(), TensorKit.QL(), TensorKit.QLpos(), TensorKit.Polar(), TensorKit.SVD(), TensorKit.SDD()) - Q, R = @constinferred leftorth(t, (3,4,2), (1,5); alg = alg) - QdQ = Q'*Q + Q, R = @constinferred leftorth(t, (3, 4, 2), (1, 5); alg=alg) + QdQ = Q' * Q @test QdQ ≈ one(QdQ) - @test Q*R ≈ permute(t, (3,4,2), (1,5)) + @test Q * R ≈ permute(t, (3, 4, 2), (1, 5)) if alg isa Polar @test isposdef(R) @test domain(R) == codomain(R) == space(t, 1)' ⊗ space(t, 5)' end end @testset "leftnull with $alg" for alg in (TensorKit.QR(), TensorKit.SVD(), TensorKit.SDD()) - N = @constinferred leftnull(t, (3,4,2),(1,5); alg = alg) - NdN = N'*N + N = @constinferred leftnull(t, (3, 4, 2), (1, 5); alg=alg) + NdN = N' * N @test NdN ≈ one(NdN) - @test norm(N'*permute(t, (3,4,2),(1,5))) < 100*eps(norm(t)) + @test norm(N' * permute(t, (3, 4, 2), (1, 5))) < 100 * eps(norm(t)) end @testset "rightorth with $alg" for alg in (TensorKit.RQ(), TensorKit.RQpos(), TensorKit.LQ(), TensorKit.LQpos(), TensorKit.Polar(), TensorKit.SVD(), TensorKit.SDD()) - L, Q = @constinferred rightorth(t, (3,4),(2,1,5); alg = alg) - QQd = Q*Q' + L, Q = @constinferred rightorth(t, (3, 4), (2, 1, 5); alg=alg) + QQd = Q * Q' @test QQd ≈ one(QQd) - @test L*Q ≈ permute(t, (3,4),(2,1,5)) + @test L * Q ≈ permute(t, (3, 4), (2, 1, 5)) if alg isa Polar @test isposdef(L) @test domain(L) == codomain(L) == space(t, 3) ⊗ space(t, 4) end end @testset "rightnull with $alg" for alg in (TensorKit.LQ(), TensorKit.SVD(), TensorKit.SDD()) - M = @constinferred rightnull(t, (3,4),(2,1,5); alg = alg) - MMd = M*M' + M = @constinferred rightnull(t, (3, 4), (2, 1, 5); alg=alg) + MMd = M * M' @test MMd ≈ one(MMd) - @test norm(permute(t, (3,4),(2,1,5))*M') < 100*eps(norm(t)) + @test norm(permute(t, (3, 4), (2, 1, 5)) * M') < 100 * eps(norm(t)) end @testset "tsvd with $alg" for alg in (TensorKit.SVD(), TensorKit.SDD()) - U, S, V = @constinferred tsvd(t, (3,4,2),(1,5); alg = alg) - UdU = U'*U + U, S, V = @constinferred tsvd(t, (3, 4, 2), (1, 5); alg=alg) + UdU = U' * U @test UdU ≈ one(UdU) - VVd = V*V' + VVd = V * V' @test VVd ≈ one(VVd) - @test U*S*V ≈ permute(t, (3,4,2),(1,5)) + @test U * S * V ≈ permute(t, (3, 4, 2), (1, 5)) end end @testset "empty tensor" begin t = TensorMap(randn, T, V1 ⊗ V2, typeof(V1)()) @testset "leftorth with $alg" for alg in (TensorKit.QR(), TensorKit.QRpos(), TensorKit.QL(), TensorKit.QLpos(), TensorKit.Polar(), TensorKit.SVD(), TensorKit.SDD()) - Q, R = @constinferred leftorth(t; alg = alg) + Q, R = @constinferred leftorth(t; alg=alg) @test Q == t @test dim(Q) == dim(R) == 0 end @testset "leftnull with $alg" for alg in (TensorKit.QR(), TensorKit.SVD(), TensorKit.SDD()) - N = @constinferred leftnull(t; alg = alg) - @test N'*N ≈ id(domain(N)) - @test N*N' ≈ id(codomain(N)) + N = @constinferred leftnull(t; alg=alg) + @test N' * N ≈ id(domain(N)) + @test N * N' ≈ id(codomain(N)) end @testset "rightorth with $alg" for alg in (TensorKit.RQ(), TensorKit.RQpos(), TensorKit.LQ(), TensorKit.LQpos(), TensorKit.Polar(), TensorKit.SVD(), TensorKit.SDD()) - L, Q = @constinferred rightorth(copy(t'); alg = alg) + L, Q = @constinferred rightorth(copy(t'); alg=alg) @test Q == t' @test dim(Q) == dim(L) == 0 end @testset "rightnull with $alg" for alg in (TensorKit.LQ(), TensorKit.SVD(), TensorKit.SDD()) - M = @constinferred rightnull(copy(t'); alg = alg) - @test M*M' ≈ id(codomain(M)) - @test M'*M ≈ id(domain(M)) + M = @constinferred rightnull(copy(t'); alg=alg) + @test M * M' ≈ id(codomain(M)) + @test M' * M ≈ id(domain(M)) end @testset "tsvd with $alg" for alg in (TensorKit.SVD(), TensorKit.SDD()) - U, S, V = @constinferred tsvd(t; alg = alg) + U, S, V = @constinferred tsvd(t; alg=alg) @test U == t @test dim(U) == dim(S) == dim(V) end @@ -410,27 +410,27 @@ for V in spacelist t = Tensor(rand, T, V1 ⊗ V1' ⊗ V2 ⊗ V2') @testset "eig and isposdef" begin - D, V = eigen(t, (1,3), (2,4)) - D̃, Ṽ = @constinferred eig(t, (1,3), (2,4)) + D, V = eigen(t, (1, 3), (2, 4)) + D̃, Ṽ = @constinferred eig(t, (1, 3), (2, 4)) @test D ≈ D̃ @test V ≈ Ṽ - VdV = V'*V - VdV = (VdV + VdV')/2 + VdV = V' * V + VdV = (VdV + VdV') / 2 @test isposdef(VdV) - t2 = permute(t, (1,3), (2,4)) - @test t2*V ≈ V*D + t2 = permute(t, (1, 3), (2, 4)) + @test t2 * V ≈ V * D @test !isposdef(t2) # unlikely for non-hermitian map - t2 = (t2 + t2'); + t2 = (t2 + t2') D, V = eigen(t2) - VdV = V'*V + VdV = V' * V @test VdV ≈ one(VdV) D̃, Ṽ = @constinferred eigh(t2) @test D ≈ D̃ @test V ≈ Ṽ - λ = minimum(minimum(real(LinearAlgebra.diag(b))) for (c,b) in blocks(D)) + λ = minimum(minimum(real(LinearAlgebra.diag(b))) for (c, b) in blocks(D)) @test isposdef(t2) == isposdef(λ) - @test isposdef(t2 - λ*one(t2) + 0.1*one(t2)) - @test !isposdef(t2 - λ*one(t2) - 0.1*one(t2)) + @test isposdef(t2 - λ * one(t2) + 0.1 * one(t2)) + @test !isposdef(t2 - λ * one(t2) - 0.1 * one(t2)) end end end @@ -439,26 +439,26 @@ for V in spacelist for p in (1, 2, 3, Inf) # Test both a normal tensor and an adjoint one. ts = (TensorMap(randn, T, V1 ⊗ V2 ⊗ V3, V4 ⊗ V5), - TensorMap(randn, T, V4 ⊗ V5, V1 ⊗ V2 ⊗ V3)') + TensorMap(randn, T, V4 ⊗ V5, V1 ⊗ V2 ⊗ V3)') for t in ts U₀, S₀, V₀, = tsvd(t) - t = rmul!(t, 1/norm(S₀, p)) - U, S, V, ϵ = @constinferred tsvd(t; trunc = truncerr(5e-1), p = p) + t = rmul!(t, 1 / norm(S₀, p)) + U, S, V, ϵ = @constinferred tsvd(t; trunc=truncerr(5e-1), p=p) # @show p, ϵ # @show domain(S) # @test min(space(S,1), space(S₀,1)) != space(S₀,1) - U′, S′, V′, ϵ′ = tsvd(t; trunc = truncerr(nextfloat(ϵ)), p = p) + U′, S′, V′, ϵ′ = tsvd(t; trunc=truncerr(nextfloat(ϵ)), p=p) @test (U, S, V, ϵ) == (U′, S′, V′, ϵ′) - U′, S′, V′, ϵ′ = tsvd(t; trunc = truncdim(ceil(Int, dim(domain(S)))), p = p) + U′, S′, V′, ϵ′ = tsvd(t; trunc=truncdim(ceil(Int, dim(domain(S)))), p=p) @test (U, S, V, ϵ) == (U′, S′, V′, ϵ′) - U′, S′, V′, ϵ′ = tsvd(t; trunc = truncspace(space(S,1)), p = p) + U′, S′, V′, ϵ′ = tsvd(t; trunc=truncspace(space(S, 1)), p=p) @test (U, S, V, ϵ) == (U′, S′, V′, ϵ′) # results with truncationcutoff cannot be compared because they don't take degeneracy into account, and thus truncate differently - U, S, V, ϵ = tsvd(t; trunc = truncbelow(1/dim(domain(S₀))), p = p) + U, S, V, ϵ = tsvd(t; trunc=truncbelow(1 / dim(domain(S₀))), p=p) # @show p, ϵ # @show domain(S) # @test min(space(S,1), space(S₀,1)) != space(S₀,1) - U′, S′, V′, ϵ′ = tsvd(t; trunc = truncspace(space(S,1)), p = p) + U′, S′, V′, ϵ′ = tsvd(t; trunc=truncspace(space(S, 1)), p=p) @test (U, S, V, ϵ) == (U′, S′, V′, ϵ′) end end @@ -471,23 +471,23 @@ for V in spacelist t = TensorMap(randn, T, W, W) s = dim(W) expt = @constinferred exp(t) - @test reshape(convert(Array, expt), (s,s)) ≈ - exp(reshape(convert(Array, t), (s,s))) + @test reshape(convert(Array, expt), (s, s)) ≈ + exp(reshape(convert(Array, t), (s, s))) @test (@constinferred sqrt(t))^2 ≈ t - @test reshape(convert(Array, sqrt(t^2)), (s,s)) ≈ - sqrt(reshape(convert(Array, t^2), (s,s))) + @test reshape(convert(Array, sqrt(t^2)), (s, s)) ≈ + sqrt(reshape(convert(Array, t^2), (s, s))) @test exp(@constinferred log(expt)) ≈ expt - @test reshape(convert(Array, log(expt)), (s,s)) ≈ - log(reshape(convert(Array, expt), (s,s))) + @test reshape(convert(Array, log(expt)), (s, s)) ≈ + log(reshape(convert(Array, expt), (s, s))) @test (@constinferred cos(t))^2 + (@constinferred sin(t))^2 ≈ id(W) - @test (@constinferred tan(t)) ≈ sin(t)/cos(t) - @test (@constinferred cot(t)) ≈ cos(t)/sin(t) + @test (@constinferred tan(t)) ≈ sin(t) / cos(t) + @test (@constinferred cot(t)) ≈ cos(t) / sin(t) @test (@constinferred cosh(t))^2 - (@constinferred sinh(t))^2 ≈ id(W) - @test (@constinferred tanh(t)) ≈ sinh(t)/cosh(t) - @test (@constinferred coth(t)) ≈ cosh(t)/sinh(t) + @test (@constinferred tanh(t)) ≈ sinh(t) / cosh(t) + @test (@constinferred coth(t)) ≈ cosh(t) / sinh(t) t1 = sin(t) @test sin(@constinferred asin(t1)) ≈ t1 @@ -512,13 +512,13 @@ for V in spacelist for T in (Float32, ComplexF64) tA = TensorMap(rand, T, V1 ⊗ V3, V1 ⊗ V3) tB = TensorMap(rand, T, V2 ⊗ V4, V2 ⊗ V4) - tA = 3//2*leftorth(tA; alg = Polar())[1] - tB = 1//5*leftorth(tB; alg = Polar())[1] + tA = 3 // 2 * leftorth(tA; alg=Polar())[1] + tB = 1 // 5 * leftorth(tB; alg=Polar())[1] tC = TensorMap(rand, T, V1 ⊗ V3, V2 ⊗ V4) t = @constinferred sylvester(tA, tB, tC) @test codomain(t) == V1 ⊗ V3 @test domain(t) == V2 ⊗ V4 - @test norm(tA*t + t*tB + tC) < (norm(tA)+norm(tB)+norm(tC))*eps(real(T))^(2/3) + @test norm(tA * t + t * tB + tC) < (norm(tA) + norm(tB) + norm(tC)) * eps(real(T))^(2 / 3) if hasfusiontensor(I) matrix(x) = reshape(convert(Array, x), dim(codomain(x)), dim(domain(x))) @test matrix(t) ≈ sylvester(matrix(tA), matrix(tB), matrix(tC)) @@ -545,8 +545,8 @@ for V in spacelist d4 = dim(domain(t2)) At = convert(Array, t) @test reshape(At, (d1, d2, d3, d4)) ≈ - reshape(convert(Array, t1), (d1, 1, d3, 1)) .* - reshape(convert(Array, t2), (1, d2, 1, d4)) + reshape(convert(Array, t1), (d1, 1, d3, 1)) .* + reshape(convert(Array, t2), (1, d2, 1, d4)) end end end @@ -555,14 +555,14 @@ for V in spacelist t1 = Tensor(rand, T, V2 ⊗ V3 ⊗ V1) t2 = Tensor(rand, T, V2 ⊗ V1 ⊗ V3) t = @constinferred (t1 ⊗ t2) - @tensor t′[1, 2, 3, 4, 5, 6] := t1[1,2,3]*t2[4,5,6] + @tensor t′[1, 2, 3, 4, 5, 6] := t1[1, 2, 3] * t2[4, 5, 6] @test t ≈ t′ end end global tf = time() printstyled("Finished tensor tests with symmetry $Istr in ", - string(round(tf-ti; sigdigits=3)), - " seconds."; bold = true, color = Base.info_color()) + string(round(tf - ti; sigdigits=3)), + " seconds."; bold=true, color=Base.info_color()) println() end @@ -580,8 +580,8 @@ end d4 = dim(domain(t2)) At = convert(Array, t) @test reshape(At, (d1, d2, d3, d4)) ≈ - reshape(convert(Array, t1), (d1, 1, d3, 1)) .* - reshape(convert(Array, t2), (1, d2, 1, d4)) + reshape(convert(Array, t1), (d1, 1, d3, 1)) .* + reshape(convert(Array, t2), (1, d2, 1, d4)) end end end