code to decompose problem for AutFN into an orbit problem
This commit is contained in:
parent
d0b1f65cad
commit
6a7b054b9f
|
|
@ -0,0 +1,306 @@
|
|||
push!(LOAD_PATH, "./")
|
||||
|
||||
using Nemo
|
||||
using Groups
|
||||
using WreathProducts
|
||||
|
||||
using GroupRings
|
||||
|
||||
import Nemo.elements
|
||||
|
||||
using JLD
|
||||
using ProgressMeter
|
||||
|
||||
function products{T<:GroupElem}(X::AbstractVector{T}, Y::AbstractVector{T})
|
||||
result = Vector{T}()
|
||||
seen = Set{T}()
|
||||
@showprogress for x in X
|
||||
for y in Y
|
||||
z = x*y
|
||||
if !in(z, seen)
|
||||
push!(seen, z)
|
||||
push!(result, z)
|
||||
end
|
||||
end
|
||||
end
|
||||
return result
|
||||
end
|
||||
|
||||
function generate_balls{T<:GroupElem}(S::Vector{T}, Id::T; radius=2)
|
||||
sizes = Vector{Int}()
|
||||
S = unshift!(S, Id)
|
||||
B = [Id]
|
||||
for i in 1:radius
|
||||
B = products(B, S);
|
||||
push!(sizes, length(B))
|
||||
end
|
||||
return B, sizes
|
||||
end
|
||||
|
||||
function elements(F::Nemo.FqNmodFiniteField)
|
||||
deg = Int(degree(F))
|
||||
char = Int(characteristic(F))
|
||||
z = (gen(F)^i for i in 0:deg-1)
|
||||
crtsn_prd = Base.product([0:char-1 for i in 1:deg]...)
|
||||
|
||||
function _it()
|
||||
for crt in crtsn_prd
|
||||
g = sum([b*a for (a,b) in zip(z,crt)])
|
||||
produce(g)
|
||||
end
|
||||
end
|
||||
return Task(_it)
|
||||
end
|
||||
|
||||
function AutFG_emb(A::AutGroup, g::WreathProductElem)
|
||||
isa(A.objectGroup, FreeGroup) || throw("Not an Aut(FN)")
|
||||
parent(g).P.n == length(A.objectGroup.gens) || throw("No natural action of $(parent(g)) on $A")
|
||||
powers = [(elt == parent(elt)() ? 0: 1) for elt in g.n.elts]
|
||||
elt = reduce(*, [A(Groups.flip_autsymbol(i))^pow for (i,pow) in enumerate(powers)])
|
||||
Groups.r_multiply!(elt, [Groups.perm_autsymbol(g.p)])
|
||||
return elt
|
||||
end
|
||||
|
||||
function (g::WreathProductElem)(a::AutGroupElem)
|
||||
g = AutFG_emb(parent(a),g)
|
||||
return g*a*g^-1
|
||||
end
|
||||
|
||||
function orbit_decomposition(G::Nemo.Group, E::Vector, rdict=GroupRings.reverse_dict(E))
|
||||
|
||||
elts = collect(elements(G))
|
||||
|
||||
tovisit = trues(E);
|
||||
orbits = Vector{Set{Int}}()
|
||||
|
||||
@showprogress "Orbit decomposition... " for i in 1:endof(E)
|
||||
if tovisit[i]
|
||||
orbit = Set{Int}()
|
||||
a = E[i]
|
||||
for g in elts
|
||||
idx = rdict[g(a)]
|
||||
tovisit[idx] = false
|
||||
push!(orbit,idx)
|
||||
end
|
||||
push!(orbits, orbit)
|
||||
end
|
||||
end
|
||||
return orbits
|
||||
end
|
||||
|
||||
function matrix_repr(g::WreathProductElem, E, E_dict)
|
||||
rep_matrix = spzeros(Int, length(E), length(E))
|
||||
for (i,e) in enumerate(E)
|
||||
j = E_dict[g(e)]
|
||||
rep_matrix[i,j] = 1
|
||||
end
|
||||
return rep_matrix
|
||||
end
|
||||
|
||||
function action_mreps(G::Nemo.Group, E, E_dict)
|
||||
result = @showprogress [matrix_repr(g, E, E_dict) for g in elements(G)]
|
||||
return result
|
||||
end
|
||||
|
||||
function chars(G::PermutationGroup)
|
||||
permtype_unsorted(σ::Nemo.perm) = [length(c) for c in cycles(σ)]
|
||||
permtype(σ::Nemo.perm) = sort(permtype_unsorted(σ))
|
||||
|
||||
χ_id(σ::Nemo.perm) = 1
|
||||
|
||||
χ_sgn(σ::Nemo.perm) = (-1)^parity(σ)
|
||||
|
||||
function χ_reg(σ::Nemo.perm)
|
||||
fixed_points = countnz([(x == y? 1 : 0) for (x,y) in enumerate(σ.d)])
|
||||
return fixed_points - 1
|
||||
end
|
||||
|
||||
χ_regsgn(σ::Nemo.perm) = (-1)^parity(σ)*χ_reg(σ)
|
||||
|
||||
function χ_regviaS3(σ::Nemo.perm)
|
||||
@assert parent(σ).n == 4
|
||||
t = permtype(σ)
|
||||
if t == [1,1,1,1]
|
||||
result = 2
|
||||
elseif t == [2,2]
|
||||
result = 2
|
||||
elseif t == [1,3]
|
||||
result = -1
|
||||
else
|
||||
result = 0
|
||||
end
|
||||
return result
|
||||
end
|
||||
|
||||
chars = [χ_id, χ_sgn, χ_regviaS3, χ_reg, χ_regsgn]
|
||||
|
||||
if G.n == 1
|
||||
return chars[1:1]
|
||||
elseif G.n == 2
|
||||
return chars[1:2]
|
||||
elseif G.n == 3
|
||||
return [chars[1:2]..., chars[4]]
|
||||
elseif G.n == 4
|
||||
return chars[1:5]
|
||||
else
|
||||
throw("Characters for $G unknown!")
|
||||
end
|
||||
end
|
||||
|
||||
function epsilon(i, g::DirectProducts.DirectProductGroupElem)
|
||||
return reduce(*, 1, ((-1)^isone(g.elts[j]) for j in 1:i))
|
||||
end
|
||||
|
||||
import Base.convert
|
||||
convert(::Type{Int}, x::Nemo.fmpz) = x.d
|
||||
|
||||
function central_projection(RG::GroupRing, char::Function, T::Type=Rational{Int})
|
||||
result = RG(T)
|
||||
for g in RG.basis
|
||||
result[g] = char(g)
|
||||
end
|
||||
return convert(T, char(RG.group())//Int(order(RG.group))*result)
|
||||
end
|
||||
|
||||
function Cstar_repr(x::GroupRingElem, matrix_reps)
|
||||
res = zeros(matrix_reps[1])
|
||||
for i in 1:length(parent(x).basis)
|
||||
res += x.coeffs[i]*matrix_reps[i]
|
||||
end
|
||||
return res
|
||||
end
|
||||
|
||||
function rankOne_projections(G::PermutationGroup, T::Type=Rational{Int})
|
||||
RG = GroupRing(G)
|
||||
projections = [central_projection(RG, χ, T) for χ in chars(G)]
|
||||
|
||||
if G.n == 1 || G.n == 2
|
||||
return projections
|
||||
elseif G.n == 3
|
||||
rankone_projs = [
|
||||
projections[1],
|
||||
projections[2],
|
||||
1//2*(one(RG) - RG(RG.group([2,1,3])))*projections[3]
|
||||
]
|
||||
return rankone_projs
|
||||
elseif G.n == 4
|
||||
rankone_projs = [
|
||||
projections[1],
|
||||
projections[2],
|
||||
1//2*(one(RG) - RG(RG.group([2,1,3,4])))*projections[3],
|
||||
1//2*(one(RG) - RG(RG.group([2,1,3,4])))*projections[4],
|
||||
1//2*(one(RG) + RG(RG.group([2,1,3,4])))*projections[5]]
|
||||
return rankone_projs
|
||||
else
|
||||
throw("Rank-one projections for $G unknown!")
|
||||
end
|
||||
|
||||
end
|
||||
|
||||
function rankOne_projections(BN::WreathProducts.WreathProduct, T::Type=Rational{Int})
|
||||
|
||||
N = BN.P.n
|
||||
# projections as elements of the group rings RSₙ
|
||||
SNprojs_nc = [rankOne_projections(PermutationGroup(i), T) for i in 1:N]
|
||||
|
||||
# embedding into group ring of BN
|
||||
RBN = GroupRing(BN)
|
||||
RFFFF_projs = [central_projection(GroupRing(BN.N), g->epsilon(i,g), T)
|
||||
for i in 0:BN.P.n]
|
||||
Qs = [RBN(q, g -> BN(g)) for q in RFFFF_projs]
|
||||
|
||||
function incl(k::Int, g::perm, WP::WreathProduct=BN)
|
||||
@assert length(g.d) + k <= WP.P.n
|
||||
arr = [1:k; g.d .+ k; (length(g.d)+k+1):WP.P.n]
|
||||
return WP(WP.P(arr))
|
||||
end
|
||||
|
||||
all_projs=[Qs[1]*RBN(p, g-> incl(0,g)) for p in SNprojs_nc[N]]
|
||||
|
||||
for i in 1:N-1
|
||||
Sk_first = [RBN(p, g->incl(0,g)) for p in SNprojs_nc[i]]
|
||||
Sk_last = [RBN(p, g->incl(i,g)) for p in SNprojs_nc[N-i]]
|
||||
append!(all_projs, [Qs[i+1]*p1*p2
|
||||
for (p1,p2) in Base.product(Sk_first,Sk_last)])
|
||||
end
|
||||
|
||||
append!(all_projs, [Qs[N+1]*RBN(p, g-> incl(0,g)) for p in SNprojs_nc[N]])
|
||||
|
||||
return all_projs
|
||||
end
|
||||
|
||||
function orthSVD(M::AbstractMatrix)
|
||||
M = full(M)
|
||||
# matrixRank = rank(M)
|
||||
fact = svdfact(M)
|
||||
sings = fact[:S]
|
||||
# @show sings[sings.>1e-14]
|
||||
M_rank = sum(fact[:S] .> maximum(size(M))*eps(eltype(fact[:S])))
|
||||
Ufactor = fact[:U]
|
||||
return Ufactor[:,1:M_rank]
|
||||
end
|
||||
|
||||
function Uπ_matrices(P_matrices; orth=orthSVD)
|
||||
U_p_matrices = Vector{Array{Float64,2}}()
|
||||
@showprogress "Computing Uπ mats..." for (i,p_mat) in enumerate(P_matrices)
|
||||
U_p = orth(p_mat)
|
||||
push!(U_p_matrices, U_p)
|
||||
end
|
||||
return U_p_matrices
|
||||
end
|
||||
|
||||
function main(N::Int)
|
||||
name = "test"
|
||||
# SOutF$(N)_E4
|
||||
|
||||
isdir(name) || mkdir(name)
|
||||
|
||||
SOutFN = AutGroup(FreeGroup(N), special=true, outer=true)
|
||||
S = generators(SOutFN);
|
||||
S = [S; [inv(s) for s in S]]
|
||||
@show S
|
||||
|
||||
info("Generating ball of radius 4")
|
||||
@time E4, sizes = generate_balls(S, SOutFN(), radius=4);
|
||||
info("Reverse dict")
|
||||
@time E_dict = GroupRings.reverse_dict(E4)
|
||||
|
||||
|
||||
BN = WreathProduct(Nemo.FiniteField(2,1, "a")[1], PermutationGroup(N))
|
||||
|
||||
if !isfile(joinpath(name, "orbits.jld"))
|
||||
info("Decomposing E into orbits of B$(N)")
|
||||
@time orbs = orbit_decomposition(BN, E4, E_dict)
|
||||
@assert sum(length(o) for o in orbs) == length(E4)
|
||||
save(joinpath(name, "orbits.jld"), "orbits", orbs)
|
||||
end
|
||||
|
||||
info("Action matrices")
|
||||
|
||||
pm = GroupRings.create_pm(E4, E_dict, sizes[2], twisted=true)
|
||||
|
||||
RSOutFN = GroupRing(SOutFN, E4, E_dict, pm)
|
||||
|
||||
PropertyT.splaplacian(RSOutFN, S)
|
||||
|
||||
E2 = E4[1:sizes[2]]
|
||||
@time BNactionE_mreps = action_mreps(BN, E2, E_dict)
|
||||
|
||||
save(joinpath(name, "matrix_reps.jld"), "matrix_reps", BNactionE_mreps)
|
||||
|
||||
info("Projections")
|
||||
@time BN_mps = rankOne_projections(BN);
|
||||
|
||||
@time π_E_projections = @showprogress [Cstar_repr(p, BNactionE_mreps) for p in BN_mps]
|
||||
|
||||
info("Uπs...")
|
||||
@time Uπs = Uπ_matrices(π_E_projections);
|
||||
|
||||
@show multiplicities = [size(U,2) for U in Uπs];
|
||||
@show dimensions = [Int(p[BN()]*Int(order(BN))) for p in BN_mps];
|
||||
@assert dot(multiplicities, dimensions) == sizes[2]
|
||||
|
||||
save(joinpath(name, "U_pis.jld"), "Uπs", Uπs, "dims", dimensions)
|
||||
end
|
||||
|
||||
@time main(4)
|
||||
Loading…
Reference in New Issue