1
0
mirror of https://github.com/kalmarek/PropertyT.jl.git synced 2024-12-26 02:30:29 +01:00

move group specific files outside of package

This commit is contained in:
kalmar 2017-03-13 15:09:03 +01:00
parent ca830c392e
commit cf38f48e6c
6 changed files with 0 additions and 718 deletions

120
AutF4.jl
View File

@ -1,120 +0,0 @@
using Combinatorics
using JuMP
import SCS: SCSSolver
import Mosek: MosekSolver
push!(LOAD_PATH, "./")
using SemiDirectProduct
using GroupAlgebras
include("property(T).jl")
const N = 4
const VERBOSE = true
function permutation_matrix(p::Vector{Int})
n = length(p)
sort(p) == collect(1:n) || throw(ArgumentError("Input array must be a permutation of 1:n"))
A = eye(n)
return A[p,:]
end
SymmetricGroup(n) = [nthperm(collect(1:n), k) for k in 1:factorial(n)]
# const SymmetricGroup = [permutation_matrix(x) for x in SymmetricGroup_perms]
function E(i, j; dim::Int=N)
@assert i≠j
k = eye(dim)
k[i,j] = 1
return k
end
function eltary_basis_vector(i; dim::Int=N)
result = zeros(dim)
if 0 < i dim
result[i] = 1
end
return result
end
v(i; dim=N) = eltary_basis_vector(i,dim=dim)
ϱ(i,j::Int,n=N) = SemiDirectProductElement(E(i,j,dim=n), v(j,dim=n))
λ(i,j::Int,n=N) = SemiDirectProductElement(E(i,j,dim=n), -v(j,dim=n))
function ɛ(i, n::Int=N)
result = eye(n)
result[i,i] = -1
return SemiDirectProductElement(result)
end
σ(permutation::Vector{Int}) =
SemiDirectProductElement(permutation_matrix(permutation))
# Standard generating set: 103 elements
function generatingset_ofAutF(n::Int=N)
indexing = [[i,j] for i in 1:n for j in 1:n if i≠j]
ϱs = [ϱ(ij...) for ij in indexing]
λs = [λ(ij...) for ij in indexing]
ɛs = [ɛ(i) for i in 1:N]
σs = [σ(perm) for perm in SymmetricGroup(n)]
S = vcat(ϱs, λs, ɛs, σs);
S = unique(vcat(S, [inv(x) for x in S]));
return S
end
#=
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!
=#
const ID = eye(N+1)
const S₁ = generatingset_ofAutF(N)
matrix_S₁ = [matrix_repr(x) for x in S₁]
const TOL=10.0^-7
matrix_S₁[1:10,:][:,1]
Δ, cm = prepare_Laplacian_and_constraints(matrix_S₁)
#solver = SCSSolver(eps=TOL, max_iters=ITERATIONS, verbose=true);
solver = MosekSolver(MSK_DPAR_INTPNT_CO_TOL_REL_GAP=TOL,
# MSK_DPAR_INTPNT_CO_TOL_PFEAS=1e-15,
# MSK_DPAR_INTPNT_CO_TOL_DFEAS=1e-15,
# MSK_IPAR_PRESOLVE_USE=0,
QUIET=!VERBOSE)
# κ, A = solve_for_property_T(S₁, solver, verbose=VERBOSE)
product_matrix = readdlm("SL₃Z.product_matrix", Int)
L = readdlm("SL₃Z.Δ.coefficients")[:, 1]
Δ = GroupAlgebraElement(L, product_matrix)
A = readdlm("matrix.A.Mosek")
κ = readdlm("kappa.Mosek")[1]
# @show eigvals(A)
@assert isapprox(eigvals(A), abs(eigvals(A)), atol=TOL)
@assert A == Symmetric(A)
const A_sqrt = real(sqrtm(A))
SOS_EOI_fp_L₁, Ω_fp_dist = check_solution(κ, A_sqrt, Δ)
κ_rational = rationalize(BigInt, κ;)
A_sqrt_rational = rationalize(BigInt, A_sqrt)
Δ_rational = rationalize(BigInt, Δ)
SOS_EOI_rat_L₁, Ω_rat_dist = check_solution(κ_rational, A_sqrt_rational, Δ_rational)

View File

@ -1,111 +0,0 @@
module AutGroups
using Groups
using Permutations
import Base: inv, ^
import Groups: IdSymbol, change_pow, GWord, ==, hash, reduce!
export AutSymbol, AutWord, GWord
export rmul_AutSymbol, lmul_AutSymbol, flip_AutSymbol, symmetric_AutSymbol
immutable AutSymbol <: GSymbol
gen::String
pow::Int
ex::Expr
end
(==)(s::AutSymbol, t::AutSymbol) = s.gen == t.gen && s.pow == t.pow
hash(s::AutSymbol, h::UInt) = hash(s.gen, hash(s.pow, hash(:AutSymbol, h)))
IdSymbol(::Type{AutSymbol}) = AutSymbol("(id)", 0, :(IdAutomorphism(N)))
function change_pow(s::AutSymbol, n::Int)
if n == 0
return one(s)
end
symbol = s.ex.args[1]
if symbol ==
return flip_AutSymbol(s.ex.args[2], pow=n)
elseif symbol == :σ
return symmetric_AutSymbol(s.ex.args[2], pow=n)
elseif symbol == :ϱ
return rmul_AutSymbol(s.ex.args[2], s.ex.args[3], pow=n)
elseif symbol ==
return lmul_AutSymbol(s.ex.args[2], s.ex.args[3], pow=n)
elseif symbol == :IdAutomorphism
return s
else
warn("Changing an unknown type of symbol! $s")
return AutSymbol(s.gen, n, s.ex)
end
end
inv(f::AutSymbol) = change_pow(f, -1*f.pow)
(^)(s::AutSymbol, n::Integer) = change_pow(s, s.pow*n)
function rmul_AutSymbol(i,j; pow::Int=1)
gen = string('ϱ',Char(8320+i), Char(8320+j)...)
return AutSymbol(gen, pow, :(ϱ($i,$j)))
end
function lmul_AutSymbol(i,j; pow::Int=1)
gen = string('λ',Char(8320+i), Char(8320+j)...)
return AutSymbol(gen, pow, :(λ($i,$j)))
end
function flip_AutSymbol(j; pow::Int=1)
gen = string('ɛ', Char(8320 + j))
return AutSymbol(gen, (2+ pow%2)%2, :(ɛ($j)))
end
function symmetric_AutSymbol(perm::Vector{Int}; pow::Int=1)
# if perm == collect(1:length(perm))
# return one(AutSymbol)
# end
perm = Permutation(perm)
ord = order(perm)
pow = pow % ord
perm = perm^pow
gen = string('σ', [Char(8320 + i) for i in array(perm)]...)
return AutSymbol(gen, 1, :(σ($(array(perm)))))
end
function getperm(s::AutSymbol)
if s.ex.args[1] == :σ
return s.ex.args[2]
else
throw(ArgumentError("$s is not a permutation automorphism!"))
end
end
typealias AutWord GWord{AutSymbol}
function simplify_perms!(W::AutWord)
reduced = true
for i in 1:length(W.symbols) - 1
current = W.symbols[i]
if current.ex.args[1] == :σ
if current.pow != 1
current = symmetric_AutSymbol(perm(current), pow=current.pow)
end
next_s = W.symbols[i+1]
if next_s.ex.args[1] == :σ
reduced = false
if next_s.pow != 1
next_s = symmetric_AutSymbol(perm(next_s), pow=next_s.pow)
end
p1 = Permutation(getperm(current))
p2 = Permutation(getperm(next_s))
W.symbols[i] = one(AutSymbol)
W.symbols[i+1] = symmetric_AutSymbol(array(p1*p2))
end
end
end
return reduced
end
end #end of module AutGroups

View File

@ -1,133 +0,0 @@
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

View File

@ -1,124 +0,0 @@
Symmetrise := function(elts)
return Unique(Concatenation(elts, List(elts, Inverse)));
end;
MYAllProducts := function(elts1, elts2)
local products, elt;
products := [];
for elt in elts1 do
products := Concatenation(products, elt*elts2);
od;
return products;
end;
Products := function(elts, n)
local products, i;
if n<=0 then
return [ ];
elif n = 1 then
return elts;
else
products := elts;
for i in [2..n] do
products := MYAllProducts(elts, products);
od;
return products;
fi;
end;
IsSupportedOn := function(basis, elt)
local elt_supp, x;
elt_supp := Support(elt);
for x in elt_supp do
if not x in basis then
return false;
fi;
od;
return true;
end;
Laplacian := function(G, generating_set)
local QG, emb, result, S, g, elt;
QG := GroupRing(Rationals, G);;
emb := Embedding(G,QG);;
S := generating_set;
result := Length(S)*One(QG);
for g in S do
result := result - g^emb;
od;
return result;
end;
Vectorise := function(elt, basis)
local result, l, i, g, coeff, axis;
Assert(0, IsSupportedOn(basis, elt),
"AssertionError: Element of interest is not supported on the basis!");
result := List(0*[1..Length(basis)]);
l := CoefficientsAndMagmaElements(elt);
for i in [1..Length(l)/2] do
g := l[2*i-1];
coeff := l[2*i];
axis := Position(basis, g);
result[axis] := result[axis] + coeff;
od;
return result;
end;
Constraints := function(basis)
local result, i, j, pos;
result := [];
for i in [1..Length(basis)] do
Add(result,[]);
od;
for i in [1..Length(basis)] do
for j in [1..Length(basis)] do
pos := Position(basis, Inverse(basis[i])*basis[j]);
if not pos = fail then
Add(result[pos], [i,j]);
fi;
od;
od;
return result;
end;
SDPGenerateAll := function(G, S, basis, name)
local QG, emb, delta, delta_sq, delta_vec, delta_sq_vec, product_constr;
Print("Initializing GroupAlgebra");
QG := GroupRing(Rationals, G);;
Print(".");
emb := Embedding(G,QG);;
Print("\n");
Print("Initializing GroupAlgebra elements: ");
delta := Laplacian(G, S);;
Print("delta! ");
delta_sq := delta^2;;
Print("delta_sq! ");
Print("\n");
Print("Check if delta_sq is supported on the given basis: ");
if not IsSupportedOn(basis, delta_sq) then
Print("delta_sq is not supported on basis\n");
return fail;
else
Print("it is!\n");
PrintTo(Concatenation("./basis.", name), basis);
Print("Written basis to ", Concatenation("./basis.", name), "\n");
delta_vec := Vectorise(delta, basis);;
PrintTo(Concatenation("./delta.", name), delta_vec);
Print("Written delta to ", Concatenation("./delta.", name), "\n");
delta_sq_vec := Vectorise(delta_sq, basis);;
PrintTo(Concatenation("./delta_sq.", name), delta_sq_vec);
Print("Written delta_sq to ", Concatenation("./delta_sq.", name), "\n");
product_constr := Constraints(basis);;
PrintTo(Concatenation("./constraints.", name), product_constr);
Print("Written Matrix Constraints to ", Concatenation("./Constraints.", name), "\n");
return "Done!";
fi;
end;;

142
SL3Z.jl
View File

@ -1,142 +0,0 @@
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)

View File

@ -1,88 +0,0 @@
module SemiDirectProduct
import Base: convert, show, isequal, ==, size, inv
import Base: +, -, *, //
export SemiDirectProductElement, matrix_repr
"""
Implements elements of a semidirect product of groups H and N, where N is normal in the product. Usually written as H N.
The multiplication inside semidirect product is defined as
(h₁, n₁) (h₂, n₂) = (h₁h₂, n₁φ(h₁)(n₂)),
where φ:H Aut(N) is a homomorphism.
In the case below we implement H = GL(n,K) and N = Kⁿ, the Affine Group (i.e. GL(n,K) Kⁿ) where elements of GL(n,K) act on vectors in Kⁿ via matrix multiplication.
# Arguments:
* `h::Array{T,2}` : square invertible matrix (element of GL(n,K))
* `n::Vector{T,1}` : vector in Kⁿ
* `φ = φ(h,n) = φ(h)(n)` :2-argument function which defines the action of GL(n,K) on Kⁿ; matrix-vector multiplication by default.
"""
immutable SemiDirectProductElement{T<:Number}
h::Array{T,2}
n::Vector{T}
φ::Function
function SemiDirectProductElement(h::Array{T,2},n::Vector{T},φ::Function)
# size(h,1) == size(h,2)|| throw(ArgumentError("h has to be square matrix"))
det(h) 0 || throw(ArgumentError("h has to be invertible!"))
new(h,n,φ)
end
end
SemiDirectProductElement{T}(h::Array{T,2}, n::Vector{T}, φ) =
SemiDirectProductElement{T}(h,n,φ)
SemiDirectProductElement{T}(h::Array{T,2}, n::Vector{T}) =
SemiDirectProductElement(h,n,*)
SemiDirectProductElement{T}(h::Array{T,2}) =
SemiDirectProductElement(h,zeros(h[:,1]))
SemiDirectProductElement{T}(n::Vector{T}) =
SemiDirectProductElement(eye(eltype(n), n))
convert{T<:Number}(::Type{T}, X::SemiDirectProductElement) =
SemiDirectProductElement(convert(Array{T,2},X.h),
convert(Vector{T},X.n),
X.φ)
size(X::SemiDirectProductElement) = (size(X.h), size(X.n))
matrix_repr{T}(X::SemiDirectProductElement{T}) =
[X.h X.n; zeros(T, 1, size(X.h,2)) [1]]
show{T}(io::IO, X::SemiDirectProductElement{T}) = print(io,
"Element of SemiDirectProduct over $T of size $(size(X)):\n",
matrix_repr(X))
function isequal{T}(X::SemiDirectProductElement{T}, Y::SemiDirectProductElement{T})
X.h == Y.h || return false
X.n == Y.n || return false
X.φ == Y.φ || return false
return true
end
function isequal{T,S}(X::SemiDirectProductElement{T}, Y::SemiDirectProductElement{S})
W = promote_type(T,S)
warn("Comparing elements with different coefficients! trying to promoting to $W")
X = convert(W, X)
Y = convert(W, Y)
return isequal(X,Y)
end
(==)(X::SemiDirectProductElement, Y::SemiDirectProductElement) = isequal(X, Y)
function semidirect_multiplication{T}(X::SemiDirectProductElement{T},
Y::SemiDirectProductElement{T})
size(X) == size(Y) || throw(ArgumentError("trying to multiply elements from different groups!"))
return SemiDirectProductElement(X.h*Y.h, X.n + X.φ(X.h, Y.n))
end
(*){T}(X::SemiDirectProductElement{T}, Y::SemiDirectProductElement{T}) =
semidirect_multiplication(X,Y)
inv{T}(X::SemiDirectProductElement{T}) =
SemiDirectProductElement(inv(X.h), X.φ(inv(X.h), -X.n))
end