Verify that topological classification is invariant under duality transformations of symmetry eigenvalues

duality-check.jl

using BandGraphs
using Crystalline
using Test
using SymmetryBases
using LinearAlgebra

@eval Crystalline begin
    struct SiteIrrep{D} <: AbstractIrrep{D}
        cdml     :: String
        g        :: SiteGroup{D}
        matrices :: Vector{Matrix{ComplexF64}}
        reality  :: Reality
        iscorep  :: Bool
        pglab    :: String
    end
    Base.position(siteir::SiteIrrep) = position(group(siteir))

    """
        siteirreps(sitegroup::SiteGroup) --> Vector{PGIrrep}

    Return the site symmetry irreps associated with the provided `SiteGroup`, derived from a
    search over isomorphic point groups.
    """
    function siteirreps(siteg::SiteGroup{D}) where D
        parent_pg, Iᵖ²ᵍ, _ = find_isomorphic_parent_pointgroup(siteg)
        pglab = label(parent_pg)
        pgirs = pgirreps(pglab, Val(D))
        
        # note that we _have to_ make a copy when re-indexing `pgir.matrices` here, since
        # .jld files apparently cache accessed content; so if we modify it, we mess with the
        # underlying data (see https://github.com/JuliaIO/JLD2.jl/issues/277)
        siteirs = map(pgirs) do pgir
            SiteIrrep{D}(label(pgir), siteg, pgir.matrices[Iᵖ²ᵍ], reality(pgir), pgir.iscorep,
                        pglab)
        end
        return siteirs
    end
    pglabel(siteir::SiteIrrep)  = siteir.pglab # associated point group label
    mulliken(siteir::SiteIrrep) = _mulliken(pglabel(siteir), label(siteir), iscorep(siteir))
end

function isproper(op)
    det(rotation(op)) > 0 ? true : false
end
siteir_name(br) = replace(br.label, "↑G"=>"")
br_name(br) = "(" * br.wyckpos * "|" * siteir_name(br) * ")"

function remap_bandreps(brs, lgirsd, wps, brs′=brs)
    @testset "SG $(brs.sgnum)" begin
        B′ = matrix(brs′, false)
        F′ = smith(B′)
        for br in brs
            n = SymVector(br, brs.irlabs, lgirsd)
            mults′ = [similar(m) for m in n.mults]
            for (i, lgirs) in enumerate(n.lgirsv)
                χs = characters(lgirs) * n.mults[i]
                χs′ = ifelse.(isproper.(group(first(lgirs))), 1, -1) .* χs
                mults′[i] = find_representation(χs′, lgirs)
            end
            nv′ = reduce(vcat, mults′)

            t = calc_detailed_topology(collect(n), B′, F′)
            t′ = calc_detailed_topology(nv′, B′, F′)
            @test t == t′
        end
    end
end

##
D = 3
timereversal = true

for sgnum in 1:MAX_SGNUM[D]
    sb, brs = compatibility_basis(sgnum, D; timereversal)
    lgirsd = lgirreps(sgnum, D)
    timereversal && (lgirsd = Dict(klab=>realify(lgirs) for (klab, lgirs) in lgirsd))
    wps = wyckoffs(sgnum, D)
    remap_bandreps(sb, lgirsd, wps, brs)
end