From 884be48c3870a2e8ded2560f578bafcd5b695c09 Mon Sep 17 00:00:00 2001 From: kalmar Date: Mon, 13 Mar 2017 16:18:42 +0100 Subject: [PATCH] Initial Commit --- AutFN.jl | 183 +++++++++++++++++++++++++++++++++++++++++++++++ GroupAlgebras.jl | 133 ++++++++++++++++++++++++++++++++++ SL3Z.jl | 142 ++++++++++++++++++++++++++++++++++++ 3 files changed, 458 insertions(+) create mode 100644 AutFN.jl create mode 100644 GroupAlgebras.jl create mode 100644 SL3Z.jl diff --git a/AutFN.jl b/AutFN.jl new file mode 100644 index 0000000..54fa1c3 --- /dev/null +++ b/AutFN.jl @@ -0,0 +1,183 @@ +using JLD +using JuMP +import SCS: SCSSolver +import Mosek: MosekSolver + +using Groups +using ProgressMeter + +#= +Note that the element + α(i,j,k) = ϱ(i,j)*ϱ(i,k)*inv(ϱ(i,j))*inv(ϱ(i,k)), +which surely belongs to ball of radius 4 in Aut(F₄) becomes trivial under the representation + Aut(F₄) → GL₄(ℤ)⋉ℤ⁴ → GL₅(ℂ). +Moreover, due to work of Potapchik and Rapinchuk [1] every real representation of Aut(Fₙ) into GLₘ(ℂ) (for m ≤ 2n-2) factors through GLₙ(ℤ)⋉ℤⁿ, so will have the same problem. + +We need a different approach: Here we actually compute in Aut(𝔽₄) +=# + +import Combinatorics.nthperm +SymmetricGroup(n) = [nthperm(collect(1:n), k) for k in 1:factorial(n)] + +function generating_set_of_AutF(N::Int) + + indexing = [[i,j] for i in 1:N for j in 1:N if i≠j] + σs = [symmetric_AutSymbol(perm) for perm in SymmetricGroup(N)[2:end]]; + ϱs = [rmul_AutSymbol(i,j) for (i,j) in indexing] + λs = [lmul_AutSymbol(i,j) for (i,j) in indexing] + ɛs = [flip_AutSymbol(i) for i in 1:N]; + + S = vcat(ϱs,λs) + S = vcat(S..., σs..., ɛs) + S = vcat(S..., [inv(g) for g in S]) + return Vector{AutWord}(unique(S)) +end + +function generating_set_of_OutF(N::Int) + indexing = [[i,j] for i in 1:N for j in 1:N if i≠j] + ϱs = [rmul_AutSymbol(i,j) for (i,j) in indexing] + λs = [lmul_AutSymbol(i,j) for (i,j) in indexing] + ɛs = [flip_AutSymbol(i) for i in 1:N]; + + S = ϱs + push!(S, λs..., ɛs...) + push!(S,[inv(g) for g in S]...) + return Vector{AutWord}(unique(S)) +end + +function generating_set_of_Sym(N::Int) + σs = [symmetric_AutSymbol(perm) for perm in SymmetricGroup(N)[2:end]]; + + S = σs + push!(S, [inv(s) for s in S]...) + return Vector{AutWord}(unique(S)) +end + + +function products(S1::Vector{AutWord}, S2::Vector{AutWord}) + result = Vector{AutWord}() + seen = Set{Vector{FGWord}}() + n = length(S1) + p = Progress(n, 1, "Computing complete products...", 50) + for (i,x) in enumerate(S1) + for y in S2 + z::AutWord = x*y + v::Vector{FGWord} = z(domain) + if !in(v, seen) + push!(seen, v) + push!(result, z) + end + end + next!(p) + end + return result +end + +function products_images(S1::Vector{AutWord}, S2::Vector{AutWord}) + result = Vector{Vector{FGWord}}() + seen = Set{Vector{FGWord}}() + n = length(S1) + + p = Progress(n, 1, "Computing images of elts in B₄...", 50) + for (i,x) in enumerate(S1) + z = x(domain) + for y in S2 + v = y(z) + if !in(v, seen) + push!(seen, v) + push!(result, v) + end + end + next!(p) + end + return result +end + +function hashed_product{T}(image::T, B, images_dict::Dict{T, Int}) + n = size(B,1) + column = zeros(Int,n) + Threads.@threads for j in 1:n + w = (B[j])(image) + k = images_dict[w] + k ≠ 0 || throw(ArgumentError( + "($i,$j): $(x^-1)*$y don't seem to be supported on basis!")) + column[j] = k + end + return column +end + +function create_product_matrix(basis::Vector{AutWord}, images) + n = length(basis) + product_matrix = zeros(Int, (n, n)); + print("Creating hashtable of images...") + @time images_dict = Dict{Vector{FGWord}, Int}(x => i + for (i,x) in enumerate(images)) + + p = Progress(n, 1, "Computing product matrix in basis...", 50) + for i in 1:n + z = (inv(basis[i]))(domain) + product_matrix[i,:] = hashed_product(z, basis, images_dict) + next!(p) + end + return product_matrix +end + +function ΔandSDPconstraints(identity::AutWord, S::Vector{AutWord}) + + println("Generating Balls of increasing radius...") + @time B₁ = vcat([identity], S) + @time B₂ = products(B₁,B₁); + @show length(B₂) + if length(B₂) != length(B₁) + @time B₃ = products(B₁, B₂) + @show length(B₃) + if length(B₃) != length(B₂) + @time B₄_images = products_images(B₁, B₃) + else + B₄_images = unique([f(domain) for f in B₃]) + end + else + B₃ = B₂ + B₄ = B₂ + B₄_images = unique([f(domain) for f in B₃]) + end + + @show length(B₄_images) + # @assert length(B₄_images) == 3425657 + + println("Creating product matrix...") + @time pm = create_product_matrix(B₂, B₄_images) + println("Creating sdp_constratints...") + @time sdp_constraints = constraints_from_pm(pm) + + L_coeff = splaplacian_coeff(S, B₂, length(B₄_images)) + Δ = GroupAlgebraElement(L_coeff, Array{Int,2}(pm)) + + return Δ, sdp_constraints +end + + +@everywhere push!(LOAD_PATH, "./") +using GroupAlgebras +include("property(T).jl") + +const symbols = [FGSymbol("x₁",1), FGSymbol("x₂",1), FGSymbol("x₃",1), FGSymbol("x₄",1), FGSymbol("x₅",1), FGSymbol("x₆",1)] + +const TOL=1e-8 +const N = 4 +const domain = Vector{FGWord}(symbols[1:N]) +const ID = one(AutWord) + +# const name = "SYM$N" +# const upper_bound=factorial(N)-TOL^(1/5) +# S() = generating_set_of_Sym(N) + +# name = "AutF$N" +# S() = generating_set_of_AutF(N) + +name = "OutF$N" +S() = generating_set_of_OutF(N) +const upper_bound=0.05 + +BLAS.set_num_threads(4) +@time check_property_T(name, ID, S; verbose=true, tol=TOL, upper_bound=upper_bound) diff --git a/GroupAlgebras.jl b/GroupAlgebras.jl new file mode 100644 index 0000000..12ebdd2 --- /dev/null +++ b/GroupAlgebras.jl @@ -0,0 +1,133 @@ +module GroupAlgebras + +import Base: convert, show, isequal, == +import Base: +, -, *, // +import Base: size, length, norm, rationalize + +export GroupAlgebraElement + + +immutable GroupAlgebraElement{T<:Number} + coefficients::AbstractVector{T} + product_matrix::Array{Int,2} + # basis::Array{Any,1} + + function GroupAlgebraElement(coefficients::AbstractVector, + product_matrix::Array{Int,2}) + + size(product_matrix, 1) == size(product_matrix, 2) || + throw(ArgumentError("Product matrix has to be square")) + new(coefficients, product_matrix) + end +end + +# GroupAlgebraElement(c,pm,b) = GroupAlgebraElement(c,pm) +GroupAlgebraElement{T}(c::AbstractVector{T},pm) = GroupAlgebraElement{T}(c,pm) + +convert{T<:Number}(::Type{T}, X::GroupAlgebraElement) = + GroupAlgebraElement(convert(AbstractVector{T}, X.coefficients), X.product_matrix) + +show{T}(io::IO, X::GroupAlgebraElement{T}) = print(io, + "Element of Group Algebra over $T of length $(length(X)):\n $(X.coefficients)") + + +function isequal{T, S}(X::GroupAlgebraElement{T}, Y::GroupAlgebraElement{S}) + if T != S + warn("Comparing elements with different coefficients Rings!") + end + X.product_matrix == Y.product_matrix || return false + X.coefficients == Y.coefficients || return false + return true +end + +(==)(X::GroupAlgebraElement, Y::GroupAlgebraElement) = isequal(X,Y) + +function add{T<:Number}(X::GroupAlgebraElement{T}, Y::GroupAlgebraElement{T}) + X.product_matrix == Y.product_matrix || throw(ArgumentError( + "Elements don't seem to belong to the same Group Algebra!")) + return GroupAlgebraElement(X.coefficients+Y.coefficients, X.product_matrix) +end + +function add{T<:Number, S<:Number}(X::GroupAlgebraElement{T}, + Y::GroupAlgebraElement{S}) + warn("Adding elements with different base rings!") + return GroupAlgebraElement(+(promote(X.coefficients, Y.coefficients)...), + X.product_matrix) +end + +(+)(X::GroupAlgebraElement, Y::GroupAlgebraElement) = add(X,Y) +(-)(X::GroupAlgebraElement) = GroupAlgebraElement(-X.coefficients, X.product_matrix) +(-)(X::GroupAlgebraElement, Y::GroupAlgebraElement) = add(X,-Y) + +function algebra_multiplication{T<:Number}(X::AbstractVector{T}, Y::AbstractVector{T}, pm::Array{Int,2}) + result = zeros(X) + for (j,y) in enumerate(Y) + if y != zero(T) + for (i, index) in enumerate(pm[:,j]) + if X[i] != zero(T) + index == 0 && throw(ArgumentError("The product don't seem to belong to the span of basis!")) + result[index] += X[i]*y + end + end + end + end + return result +end + +function group_star_multiplication{T<:Number}(X::GroupAlgebraElement{T}, + Y::GroupAlgebraElement{T}) + X.product_matrix == Y.product_matrix || ArgumentError( + "Elements don't seem to belong to the same Group Algebra!") + result = algebra_multiplication(X.coefficients, Y.coefficients, X.product_matrix) + return GroupAlgebraElement(result, X.product_matrix) +end + +function group_star_multiplication{T<:Number, S<:Number}( + X::GroupAlgebraElement{T}, + Y::GroupAlgebraElement{S}) + S == T || warn("Multiplying elements with different base rings!") + return group_star_multiplication(promote(X,Y)...) +end + +(*){T<:Number, S<:Number}(X::GroupAlgebraElement{T}, + Y::GroupAlgebraElement{S}) = group_star_multiplication(X,Y); + +(*){T<:Number}(a::T, X::GroupAlgebraElement{T}) = GroupAlgebraElement( + a*X.coefficients, X.product_matrix) + +function scalar_multiplication{T<:Number, S<:Number}(a::T, + X::GroupAlgebraElement{S}) + promote_type(T,S) == S || warn("Scalar and coefficients are in different rings! Promoting result to $(promote_type(T,S))") + return GroupAlgebraElement(a*X.coefficients, X.product_matrix) +end + +(*){T<:Number}(a::T,X::GroupAlgebraElement) = scalar_multiplication(a, X) + +//{T<:Rational, S<:Rational}(X::GroupAlgebraElement{T}, a::S) = + GroupAlgebraElement(X.coefficients//a, X.product_matrix) + +//{T<:Rational, S<:Integer}(X::GroupAlgebraElement{T}, a::S) = + X//convert(T,a) + +length(X::GroupAlgebraElement) = length(X.coefficients) +size(X::GroupAlgebraElement) = size(X.coefficients) + +function norm(X::GroupAlgebraElement, p=2) + if p == 1 + return sum(abs(X.coefficients)) + elseif p == Inf + return max(abs(X.coefficients)) + else + return norm(X.coefficients, p) + end +end + +ɛ(X::GroupAlgebraElement) = sum(X.coefficients) + +function rationalize{T<:Integer, S<:Number}( + ::Type{T}, X::GroupAlgebraElement{S}; tol=eps(S)) + v = rationalize(T, X.coefficients, tol=tol) + return GroupAlgebraElement(v, X.product_matrix) +end + +end diff --git a/SL3Z.jl b/SL3Z.jl new file mode 100644 index 0000000..68a867f --- /dev/null +++ b/SL3Z.jl @@ -0,0 +1,142 @@ +using JLD +using JuMP +import Primes: isprime +import SCS: SCSSolver +import Mosek: MosekSolver + +using Mods + +using Groups + +function SL_generatingset(n::Int) + + indexing = [(i,j) for i in 1:n for j in 1:n if i≠j] + + S = [E(i,j,N=n) for (i,j) in indexing]; + S = vcat(S, [convert(Array{Int,2},x') for x in S]); + S = vcat(S, [convert(Array{Int,2},inv(x)) for x in S]); + return unique(S) +end + +function E(i::Int, j::Int; val=1, N::Int=3, mod=Inf) + @assert i≠j + m = eye(Int, N) + m[i,j] = val + if mod == Inf + return m + else + return [Mod(x,mod) for x in m] + end +end + +function cofactor(i,j,M) + z1 = ones(Bool,size(M,1)) + z1[i] = false + + z2 = ones(Bool,size(M,2)) + z2[j] = false + + return M[z1,z2] +end + +import Base.LinAlg.det + +function det(M::Array{Mod,2}) + if size(M,1) ≠ size(M,2) + d = Mod(0,M[1,1].mod) + elseif size(M,1) == 2 + d = M[1,1]*M[2,2] - M[1,2]*M[2,1] + else + d = zero(eltype(M)) + for i in 1:size(M,1) + d += (-1)^(i+1)*M[i,1]*det(cofactor(i,1,M)) + end + end +# @show (M, d) + return d +end + +function adjugate(M) + K = similar(M) + for i in 1:size(M,1), j in 1:size(M,2) + K[j,i] = (-1)^(i+j)*det(cofactor(i,j,M)) + end + return K +end + +import Base: inv, one, zero, * + +one(::Type{Mod}) = 1 +zero(::Type{Mod}) = 0 +zero(x::Mod) = Mod(x.mod) + +function inv(M::Array{Mod,2}) + d = det(M) + d ≠ 0*d || thow(ArgumentError("Matrix is not invertible!")) + return inv(det(M))*adjugate(M) + return adjugate(M) +end + +function SL_generatingset(n::Int, p::Int) + (p > 1 && n > 1) || throw(ArgumentError("Both n and p should be integers!")) + isprime(p) || throw(ArgumentError("p should be a prime number!")) + + indexing = [(i,j) for i in 1:n for j in 1:n if i≠j] + S = [E(i,j, N=n, mod=p) for (i,j) in indexing] + S = vcat(S, [inv(s) for s in S]) + S = vcat(S, [permutedims(x, [2,1]) for x in S]); + + return unique(S) +end + +function products{T}(U::AbstractVector{T}, V::AbstractVector{T}) + result = Vector{T}() + for u in U + for v in V + push!(result, u*v) + end + end + return unique(result) +end + +function ΔandSDPconstraints(identity, S) + B₁ = vcat([identity], S) + B₂ = products(B₁, B₁); + B₃ = products(B₁, B₂); + B₄ = products(B₁, B₃); + @assert B₄[1:length(B₂)] == B₂ + + product_matrix = create_product_matrix(B₄,length(B₂)); + sdp_constraints = constraints_from_pm(product_matrix, length(B₄)) + L_coeff = splaplacian_coeff(S, B₂, length(B₄)); + Δ = GroupAlgebraElement(L_coeff, product_matrix) + + return Δ, sdp_constraints +end + + + + + + +@everywhere push!(LOAD_PATH, "./") +using GroupAlgebras +include("property(T).jl") + +const N = 3 + +const name = "SL$(N)Z" +const ID = eye(Int, N) +S() = SL_generatingset(N) +const upper_bound=0.27 + + +# const p = 7 +# const upper_bound=0.738 # (N,p) = (3,7) + +# const name = "SL($N,$p)" +# const ID = [Mod(x,p) for x in eye(Int,N)] +# S() = SL_generatingset(N, p) + +BLAS.set_num_threads(4) +@time check_property_T(name, ID, S; verbose=true, tol=1e-10, upper_bound=upper_bound)