From 599b6e3b2646f3acaffc5d999ac31cefa8c654c7 Mon Sep 17 00:00:00 2001 From: Thomas Christensen Date: Fri, 13 Sep 2024 15:54:54 +0200 Subject: [PATCH] steps toward a new `BandRep` struct (called `NewBandRep` for now) - this also removes some fields from the original `BandRep`, namely the `decomposable` and `allpaths` fields; they are not used anywhere and cannot be easily recomputed in `calc_bandrep`. --- build/write_bandreps_from_bandrepset.jl | 6 - src/Crystalline.jl | 5 +- src/bandrep.jl | 21 ++-- src/calc_bandreps.jl | 64 +++++----- src/show.jl | 133 +++++--------------- src/types.jl | 10 +- src/types_symmetry_vectors.jl | 158 ++++++++++++++++++++++++ test/calc_bandreps.jl | 6 +- 8 files changed, 245 insertions(+), 158 deletions(-) create mode 100644 src/types_symmetry_vectors.jl diff --git a/build/write_bandreps_from_bandrepset.jl b/build/write_bandreps_from_bandrepset.jl index 5649e616..4364577a 100644 --- a/build/write_bandreps_from_bandrepset.jl +++ b/build/write_bandreps_from_bandrepset.jl @@ -28,12 +28,6 @@ function bandrep2csv(io::IO, BRS::BandRepSet) end println(io) - print(io, "Decomposable") - for br in BRS - print(io, "|", br.decomposable) - end - println(io) - irlabs = BRS.irlabs for (idx, (klab,kv)) in enumerate(zip(BRS.klabs, BRS.kvs)) print(io, klab, ":") diff --git a/src/Crystalline.jl b/src/Crystalline.jl index 8d3595e2..c0e4cbe9 100644 --- a/src/Crystalline.jl +++ b/src/Crystalline.jl @@ -73,6 +73,9 @@ export SymOperation, # types irreplabels, klabels, # ::BandRep & ::BandRepSet isspinful +include("types_symmetry_vectors.jl") +export SymmetryVector, NewBandRep, Collection + include("notation.jl") export schoenflies, iuc, centering, seitz, mulliken @@ -103,7 +106,7 @@ include("assembly/generators/spacegroup.jl") include("assembly/generators/subperiodicgroup.jl") export generate, generators -include("show.jl") # custom printing for structs defined in src/types.jl +include("show.jl") # printing of structs from src/[types.jl, types_symmetry_vectors.jl] include("orders.jl") diff --git a/src/bandrep.jl b/src/bandrep.jl index da6c5f0d..92719544 100644 --- a/src/bandrep.jl +++ b/src/bandrep.jl @@ -20,10 +20,15 @@ function array2struct(M::Matrix{String}, sgnum::Integer, allpaths::Bool=false, s temp .= split_paren.(@view M[2,2:end]) # same size, so reuse array label, dim = getindex.(temp, 1), parse.(Int, getindex.(temp, 2)) # label of bandrep - - decomposable = parse.(Bool, vec(@view M[3,2:end])) # whether bandrep can be further decomposed - - brtags = collect(eachcol(@view M[4:end, 2:end])) # set of irreps that jointly make up the bandrep + + # whether M contains info on decomposability; we don't use this anymore, but need to + # know to parse the rest of the contents correctly (we used to always include this info + # but might not in the future; so protect against this) + has_decomposable_info = M[3,1] == "Decomposable" + # decomposable = parse.(Bool, vec(@view M[3,2:end])) # whether BR can be BR-decomposed + + # set of irreps that jointly make up the bandrep + brtags = collect(eachcol(@view M[3+has_decomposable_info:end, 2:end])) for br in brtags br .= replace.(br, Ref(r"\([1-9]\)"=>"")) # get rid of irrep dimension info end @@ -34,15 +39,15 @@ function array2struct(M::Matrix{String}, sgnum::Integer, allpaths::Bool=false, s else # single-valued irreps only (spinless systems) delidxs = findall(map(isspinful, brtags)) end - for vars in (brtags, wyckpos, sitesym, label, dim, decomposable) + for vars in (brtags, wyckpos, sitesym, label, dim) deleteat!(vars, delidxs) end irlabs, irvecs = get_irrepvecs(brtags) - BRs = BandRep.(wyckpos, sitesym, label, dim, decomposable, map(isspinful, brtags), - irvecs, Ref(irlabs)) + BRs = BandRep.(wyckpos, sitesym, label, dim, map(isspinful, brtags), irvecs, + Ref(irlabs)) - return BandRepSet(sgnum, BRs, kvs, klabs, irlabs, allpaths, spinful, timereversal) + return BandRepSet(sgnum, BRs, kvs, klabs, irlabs, spinful, timereversal) end diff --git a/src/calc_bandreps.jl b/src/calc_bandreps.jl index 52c9e2fd..c1227d0b 100644 --- a/src/calc_bandreps.jl +++ b/src/calc_bandreps.jl @@ -98,7 +98,6 @@ Return the band representation induced by the provided `SiteIrrep` evaluated at a `SymOperation` `h`. """ function induce_bandrep(siteir::SiteIrrep{D}, h::SymOperation{D}, kv::KVec{D}) where D - siteg = group(siteir) wp = position(siteg) kv′ = constant(h*kv) # <-- TODO: Why only constant part? @@ -151,27 +150,35 @@ function subduce_onto_lgirreps( return m′ end +# ---------------------------------------------------------------------------------------- # + function calc_bandrep( - siteir::SiteIrrep{D}, - lgirsd::Dict{String, <:AbstractVector{LGIrrep{D}}}; - irlabs::Vector{String} = reduce_dict_of_vectors(label, lgirsd), - irdims::Vector{Int} = reduce_dict_of_vectors(irdim, lgirsd) + siteir :: SiteIrrep{D}, + lgirsv :: AbstractVector{<:AbstractVector{LGIrrep{D}}}, + timereversal :: Bool ) where D - irvec = Int[] - for (_, lgirs) in lgirsd - kv_array = subduce_onto_lgirreps(siteir, lgirs) - append!(irvec, Int.(kv_array)) - end + multsv = [subduce_onto_lgirreps(siteir, lgirs) for lgirs in lgirsv] - irdims_Γmask = [irlab[1]=='Γ' for irlab in irlabs] .* irdims # masked w/ Γ-irreps only - brdim = dot(irvec, irdims_Γmask) - wycklab = label(position(group(siteir))) - spinful = false # NB: default; Crystalline currently doesn't have spinful irreps - decomposable = false # NB: placeholder, because we don't know the true answer presently - - return BandRep(wycklab, siteir.pglabel, mulliken(siteir)*"↑G", brdim, decomposable, - spinful, irvec, irlabs) + occupation = sum(zip(first(multsv), first(lgirsv)); init=0) do (m, lgir) + m * irdim(lgir) + end + n = SymmetryVector(lgirsv, multsv, occupation) + + spinful = false # NB: default; Crystalline currently doesn't have spinful irreps + + return NewBandRep(siteir, n, timereversal, spinful) +end +function calc_bandrep( + siteir :: SiteIrrep{D}; + timereversal :: Bool=true, + allpaths :: Bool=false + ) where D + lgirsd = lgirreps(num(siteir), Val(D)) + allpaths || filter!(((_, lgirs),) -> isspecial(first(lgirs)), lgirsd) + timereversal && realify!(lgirsd) + lgirsv = [lgirs for lgirs in values(lgirsd)] + return calc_bandrep(siteir, lgirsv, timereversal) end # ---------------------------------------------------------------------------------------- # @@ -212,29 +219,20 @@ function calc_bandreps( # get all the little group irreps that we want to subduce onto lgirsd = lgirreps(sgnum, Val(D)) - allpaths || filter!(((_, lgirs),) -> isspecial(first(lgirs)), lgirsd) + allpaths || filter!(((_, lgirs),) -> isspecial(first(lgirs)), lgirsd) timereversal && realify!(lgirsd) - - irlabs = reduce_dict_of_vectors(label, lgirsd) - irdims = reduce_dict_of_vectors(irdim, lgirsd) + lgirsv = [lgirs for lgirs in values(lgirsd)] # get the bandreps induced by every maximal site symmetry irrep wps = wyckoffs(sgnum, Dᵛ) sg = spacegroup(sgnum, Dᵛ) sitegs = findmaximal(sitegroup.(Ref(sg), wps)) - brs = BandRep[] + brs = NewBandRep{D}[] for siteg in sitegs - siteirs = siteirreps(siteg) + siteirs = siteirreps(siteg; mulliken=true) timereversal && (siteirs = realify(siteirs)) - append!(brs, calc_bandrep.(siteirs, Ref(lgirsd); irlabs=irlabs, irdims=irdims)) + append!(brs, calc_bandrep.(siteirs, Ref(lgirsv), Ref(timereversal))) end - # TODO: There's potentially some kind of bad implicit dependence/overreliance on fixed - # on iteration order of `Dict`s here and related methods; probably OK, but unwise - klabs = collect(keys(lgirsd)) - kvs = position.(first.(values(lgirsd))) - spinful = false # NB: default; Crystalline currently doesn't have spinful irreps - decomposable = false # NB: placeholder, because we don't know the true answer presently - - return BandRepSet(sgnum, brs, kvs, klabs, irlabs, allpaths, spinful, timereversal) + return Collection(brs) end \ No newline at end of file diff --git a/src/show.jl b/src/show.jl index d7472947..38ecee58 100644 --- a/src/show.jl +++ b/src/show.jl @@ -527,109 +527,42 @@ function show(io::IO, ::MIME"text/plain", BRS::BandRepSet) ) # print k-vec labels - print(io, " KVecs (", hasnonmax(BRS) ? "incl. non-maximal" : "maximal only", "): ") + print(io, " KVecs: ") join(io, klabels(BRS), ", ") - - # EARLIER MANUAL LAYOUTS: DISCONTINUED - - #= - # === EBRS-by-column layout === - # prep-work to figure out how many bandreps we can write to the io - cols_avail = displaysize(io)[2] # available cols in io (cannot write to all of it; subtract 2) - indent = 3 - μ_maxdigs = maximum(ndigits∘dim, BRS) # implicitly assuming this to be ≤2... - maxcols_irr = maximum(length, irreplabels(BRS)) - cols_brlabs = length.(label.(BRS)) - cumcols_brlabs = cumsum(cols_brlabs .+ 1) - cols_requi = cumcols_brlabs .+ (indent + maxcols_irr + 2) - @show cols_requi - room_for = findlast(≤(cols_avail), cols_requi) - abbreviate = room_for < length(BRS) ? true : false - rowend = abbreviate ? '…' : '║' - - # print EBR header - print(io, ' '^(indent+maxcols_irr+1), '║') - for idxᵇʳ in 1:room_for - br = BRS[idxᵇʳ] - print(io, ' ', chop(label(br), tail=2),' ', idxᵇʳ ≠ room_for ? '│' : rowend) - end - println(io) +end - # print irrep content, line by line - for idxⁱʳʳ in 1:Nⁱʳʳ - irlab = irreplabels(BRS)[idxⁱʳʳ] - print(io, ' '^indent, irlab, ' '^(maxcols_irr + 1 - length(irlab)), '║') - for idxᵇʳ in 1:room_for - x = BRS[idxᵇʳ][idxⁱʳʳ] - addspace = div(cols_brlabs[idxᵇʳ], 2) - print(io, ' '^addspace) - iszero(x) ? print(io, '·') : print(io, x) - print(io, ' '^(cols_brlabs[idxᵇʳ] - addspace-1)) - print(io, idxᵇʳ ≠ room_for ? '│' : rowend) - end - println(io) - end +# ---------------------------------------------------------------------------------------- # +# SymmetryVector - # print band filling - print(io, ' '^indent, 'μ', ' '^maxcols_irr, '║') - for idxᵇʳ in 1:room_for - μ = dim(BRS[idxᵇʳ]) - addspace = div(cols_brlabs[idxᵇʳ], 2)+1 - print(io, ' '^(addspace-ndigits(μ)), μ) - print(io, ' '^(cols_brlabs[idxᵇʳ] - addspace)) - print(io, idxᵇʳ ≠ room_for ? '│' : rowend) - end - =# - - #= - # === EBRs-by-rows layout === - μ_maxdigs = maximum(ndigits∘dim, BRS) - cols_brlab = maximum(x->length(label(x)), BRS)+1 - cols_irstart = cols_brlab+4 - cols_avail = displaysize(io)[2]-2 # available cols in io (cannot write to all of it; subtract 2) - cols_requi = sum(x->length(x)+3, irreplabels(BRS))+cols_irstart+μ_maxdigs+3 # required cols for irrep labels & band reps - if cols_requi > cols_avail - cols_toomany = ceil(Int, (cols_requi-cols_avail)/2) + 2 # +2 is to make room for ' … ' extender - cols_midpoint = div(cols_requi-cols_irstart,2)+cols_irstart - cols_skipmin = cols_midpoint - cols_toomany - cols_skipmax = cols_midpoint + cols_toomany - cols_eachstart = [0; cumsum(length.(irreplabels(BRS)).+3)].+cols_irstart - iridx_skiprange = [idx for (idx, col_pos) in enumerate(cols_eachstart) if cols_skipmin ≤ col_pos ≤ cols_skipmax] - abbreviate = true - else - abbreviate = false +function Base.show(io :: IO, ::MIME"text/plain", n :: SymmetryVector) + print(io, length(n)-1, "-irrep ", typeof(n), ":\n ") + show(io, n) +end +function Base.show(io :: IO, n :: SymmetryVector) + print(io, "[") + for (i, (mults_k, lgirs_k)) in enumerate(zip(n.multsv, n.lgirsv)) + printstyled(io, + Crystalline.symvec2string(mults_k, label.(lgirs_k); braces=false); + color=iseven(i) ? :normal : :light_blue) + i ≠ length(n.multsv) && print(io, ", ") end + print(io, "]") + printstyled(io, " (", n.occupation, " band", n.occupation ≠ 1 ? "s" : "", ")"; + color=:light_black) +end - print(io, " "^(cols_irstart-1),'║'); # align with spaces - for (iridx,lab) in enumerate(irreplabels(BRS)) # irrep labels - if abbreviate && iridx ∈ iridx_skiprange - if iridx == first(iridx_skiprange) - print(io, "\b … ") - end - else - print(io, ' ', lab, " │") - end - end - println(io, ' '^μ_maxdigs, "μ", " ║") # band-filling column header - # print each bandrep - for (bridx,BR) in enumerate(BRS) - μ = dim(BR) - print(io, " ", label(BR), # bandrep label - " "^(cols_brlab-length(label(BR))), '║') - for (iridx,x) in enumerate(BR) # iterate over vector representation of band rep - if abbreviate && iridx ∈ iridx_skiprange - if iridx == first(iridx_skiprange) - print(io, mod(bridx,4) == 0 ? "\b … " : "\b ") - end - else - print(io, " ") - !iszero(x) ? print(io, x) : print(io, '·') - print(io, " "^(length(irreplabels(BRS)[iridx])-1), '│') # assumes we will never have ndigit(x) != 1 - end - end - - print(io, ' '^(1+μ_maxdigs-ndigits(μ)), μ, " ║") # band-filling - if bridx != length(BRS); println(io); end - end - =# +# ---------------------------------------------------------------------------------------- # +# NewBandRep + +function Base.show(io :: IO, ::MIME"text/plain", br :: NewBandRep) + print(io, length(br.n)-1, "-irrep ", typeof(br), ":\n ") + print(io, "(", ) + printstyled(io, label(position(br.siteir)); bold=true) + print(io, "|") + printstyled(io, label(br.siteir); bold=true) + print(io, "): ") + show(io, br.n) +end +function Base.show(io :: IO, br :: NewBandRep) + print(io, "(", label(position(br.siteir)), "|", label(br.siteir), ")") end \ No newline at end of file diff --git a/src/types.jl b/src/types.jl index 5f477d12..75e13b60 100644 --- a/src/types.jl +++ b/src/types.jl @@ -906,11 +906,10 @@ end $(TYPEDEF)$(TYPEDFIELDS) """ struct BandRep <: AbstractVector{Int} - wyckpos::String # Wyckoff position that induces the BR (TODO: type as `::WyckoffPosition` instead of `::String`) + wyckpos::String # Wyckoff position that induces the BR sitesym::String # Site-symmetry point group of Wyckoff pos (IUC notation) label::String # Symbol ρ↑G, with ρ denoting the irrep of the site-symmetry group dim::Int # Dimension (i.e. # of bands) in band rep - decomposable::Bool # Whether a given bandrep can be decomposed further spinful::Bool # Whether a given bandrep involves spinful irreps ("\bar"'ed irreps) irvec::Vector{Int} # Vector that references irlabs of a parent BandRepSet; nonzero # entries correspond to an element in the band representation @@ -925,11 +924,8 @@ irreplabels(BR::BandRep) = BR.irlabs dim(BR::BandRep) --> Int Return the number of bands included in the provided `BandRep`. - -If the bands are "nondetachable" (i.e. if `BR.decomposable = false`), this is equal to a -band connectivity μ. """ -dim(BR::BandRep) = BR.dim +dim(BR::BandRep) = BR.dim # TODO: Deprecate to `occupation` instead # define the AbstractArray interface for BandRep size(BR::BandRep) = (size(BR.irvec)[1] + 1,) # number of irreps sampled by BandRep + 1 (filling) @@ -956,13 +952,11 @@ struct BandRepSet <: AbstractVector{BandRep} kvs::Vector{<:KVec} # Vector of 𝐤-points # TODO: Make parametric klabs::Vector{String} # Vector of associated 𝐤-labels (in CDML notation) irlabs::Vector{String} # Vector of (sorted) CDML irrep labels at _all_ 𝐤-points - allpaths::Bool # Whether all paths (true) or only maximal 𝐤-points (false) are included spinful::Bool # Whether the band rep set includes (true) or excludes (false) spinful irreps timereversal::Bool # Whether the band rep set assumes time-reversal symmetry (true) or not (false) end num(BRS::BandRepSet) = BRS.sgnum klabels(BRS::BandRepSet) = BRS.klabs -hasnonmax(BRS::BandRepSet) = BRS.allpaths irreplabels(BRS::BandRepSet) = BRS.irlabs isspinful(BRS::BandRepSet) = BRS.spinful reps(BRS::BandRepSet) = BRS.bandreps diff --git a/src/types_symmetry_vectors.jl b/src/types_symmetry_vectors.jl new file mode 100644 index 00000000..db8d7ced --- /dev/null +++ b/src/types_symmetry_vectors.jl @@ -0,0 +1,158 @@ +# ---------------------------------------------------------------------------------------- # +# AbstractSymmetryVector +abstract type AbstractSymmetryVector{D} <: AbstractVector{Int} end + +# ::: API ::: +# multiplicities(::AbstractSymmetryVector) --> AbstractVector{<:AbstractVector{Int}} +# irreps(::AbstractSymmetryVector) --> AbstractVector{<:IrrepCollection} +# occupation(::AbstractSymmetryVector) --> Int + +# ::: AbstractArray interface ::: +Base.size(n::AbstractSymmetryVector) = (mapreduce(length, +, multiplicities(n)) + 1,) +@propagate_inbounds function Base.getindex(n::AbstractSymmetryVector, i::Int) + Nⁱʳ = length(n) + @boundscheck i > Nⁱʳ && throw(BoundsError(n, i)) + i == Nⁱʳ && return occupation(n) + i₀ = i₁ = 0 + for mults in multiplicities(n) + i₁ += length(mults) + if i ≤ i₁ + return mults[i-i₀] + end + i₀ = i₁ + end + error("unreachable reached") +end +Base.iterate(n::AbstractSymmetryVector) = multiplicities(n)[1][1], (1, 1) +function Base.iterate(n::AbstractSymmetryVector, state::Tuple{Int, Int}) + i, j = state + j += 1 + if i > length(multiplicities(n)) + return nothing + end + if j > length(multiplicities(n)[i]) + if i == length(multiplicities(n)) + return occupation(n), (i+1, j) + end + i += 1 + j = 1 + end + return multiplicities(n)[i][j], (i, j) +end + +# ---------------------------------------------------------------------------------------- # +# SymmetryVector + +mutable struct SymmetryVector{D} <: AbstractSymmetryVector{D} + const lgirsv :: Vector{IrrepCollection{LGIrrep{D}}} + const multsv :: Vector{Vector{Int}} + occupation :: Int +end +# ::: AbstractSymmetryVector interface ::: +irreps(n::SymmetryVector) = n.lgirsv +multiplicities(n::SymmetryVector) = n.multsv +occupation(n::SymmetryVector) = n.occupation + +# ::: AbstractArray interface beyond AbstractSymmetryVector ::: +function Base.similar(n::SymmetryVector{D}) where D + SymmetryVector{D}(n.lgirsv, similar(n.multsv), 0) +end +@propagate_inbounds function Base.setindex!(n::SymmetryVector, v::Int, i::Int) + Nⁱʳ = length(n) + @boundscheck i > Nⁱʳ && throw(BoundsError(n, i)) + i == Nⁱʳ && (n.occupation = v; return v) + i₀ = i₁ = 0 + for mults in multiplicities(n) + i₁ += length(mults) + if i ≤ i₁ + mults[i-i₀] = v + return v + end + i₀ = i₁ + end +end + +# ::: Optimizations and utilities ::: +function Base.Vector(n::SymmetryVector) + nv = Vector{Int}(undef, length(n)) + i = 1 + for mults in n.multsv + N = length(mults) + copyto!(nv, i, mults, 1, N) + i += N + end + nv[end] = n.occupation + return nv +end +irreplabels(n::SymmetryVector) = [label(ir) for ir in Iterators.flatten(n.lgirsv)] +klabels(n::SymmetryVector) = [klabel(first(irs)) for irs in n.lgirsv] + +# ---------------------------------------------------------------------------------------- # +# NewBandRep + +struct NewBandRep{D} <: AbstractSymmetryVector{D} + siteir :: SiteIrrep{D} + n :: SymmetryVector{D} + timereversal :: Bool + spinful :: Bool +end + +# ::: AbstractSymmetryVector interface ::: +irreps(br::NewBandRep) = irreps(br.n) +multiplicities(br::NewBandRep) = multiplicities(br.n) +occupation(br::NewBandRep) = occupation(br.n) + +# ::: AbstractArray interface beyond AbstractSymmetryVector ::: +Base.setindex!(br::NewBandRep{D}, v::Int, i::Int) where D = (br.n[i] = v) +function Base.similar(br::NewBandRep{D}) where D + NewBandRep{D}(br.siteir, similar(br.n), br.timereversal, br.spinful) +end +Base.Vector(br::NewBandRep) = Vector(br.n) + +# ::: Utilities ::: +num(br::NewBandRep) = num(br.siteir) +irreplabels(br::NewBandRep) = irreplabels(br.n) +klabels(br::NewBandRep) = klabels(br.n) + +# ::: Conversion to BandRep ::: +function Base.convert(::Type{BandRep}, br::NewBandRep{D}) where D + wyckpos = label(position(br.siteir)) + sitesym = br.siteir.pglabel + siteirlabel = label(br.siteir)*"↑G" + dim = occupation(br) + spinful = br.spinful + irvec = collect(br)[1:end-1] + irlabs = irreplabels(br) + return BandRep(wyckpos, sitesym, siteirlabel, dim, spinful, irvec, irlabs) +end + +# ---------------------------------------------------------------------------------------- # +# Collection +struct Collection{T} <: AbstractVector{T} + vs :: Vector{T} +end + +# ::: AbstractArray interface ::: +Base.size(c::Collection) = size(c.vs) +Base.getindex(c::Collection, i::Int) = c.vs[i] +Base.setindex!(c::Collection{T}, v::T, i::Int) where T = (c.vs[i]=v) +Base.similar(c::Collection{T}) where T = Collection{T}(similar(c.vs)) +Base.iterate(c::Collection) = iterate(c.vs) +Base.iterate(c::Collection, state) = iterate(c.vs, state) + +# ::: Utilities ::: +num(brs::Collection{NewBandRep}) = num(first(brs).siteir) +irreplabels(brs::Collection{NewBandRep}) = irreplabels(first(brs).n) +klabels(brs::Collection{NewBandRep}) = klabels(first(brs).n) + +# ::: Conversion to BandRepSet ::: +function Base.convert(::Type{BandRepSet}, brs::Collection{NewBandRep{D}}) where D + sgnum = num(first(brs)) + bandreps = convert.(Ref(BandRep), brs) + kvs = [position(first(lgirs)) for lgirs in irreps(first(brs))] + klabs = klabels(first(brs)) + irlabs = irreplabels(first(brs)) + spinful = first(brs).spinful + timereversal = first(brs).timereversal + return BandRepSet(sgnum, bandreps, kvs, klabs, irlabs, spinful, timereversal) +end \ No newline at end of file diff --git a/test/calc_bandreps.jl b/test/calc_bandreps.jl index e6c43bc3..c4fc0f50 100644 --- a/test/calc_bandreps.jl +++ b/test/calc_bandreps.jl @@ -16,7 +16,8 @@ using Crystalline: dlm2struct # by the latter; basically check that the data behind the former is up-to-date) for sgnum in 1:17 for timereversal in (false, true) - @test calc_bandreps(sgnum, Val(2); timereversal, allpaths = false) == + @test convert.(BandRep, + calc_bandreps(sgnum, Val(2); timereversal, allpaths = false)) == bandreps(sgnum, 2; timereversal, allpaths = false) end end @@ -30,7 +31,8 @@ end for timereversal in (false, true) had_tr_error = false - brsᶜ = calc_bandreps(sgnum, Val(3); timereversal, allpaths = false) + _brsᶜ = calc_bandreps(sgnum, Val(3); timereversal, allpaths = false) + brsᶜ = convert(BandRepSet, _brsᶜ) brsʳ = bandreps(sgnum, 3; timereversal, allpaths = false) # find a permutation of the irreps in `brsʳ` that matches `brsᶜ`'s sorting, s.t.