Initial (working) code

This commit is contained in:
kalmar 2016-12-19 15:44:52 +01:00
commit 39f0b86af2
3 changed files with 301 additions and 0 deletions

114
GroupAlgebras.jl Normal file
View File

@ -0,0 +1,114 @@
module GroupAlgebras
import Base: convert, show, isequal, ==
import Base: +, -, *, //
import Base: size, length, norm
export GroupAlgebraElement
immutable GroupAlgebraElement{T<:Number}
coordinates::Vector{T}
product_matrix::Array{Int,2}
# basis::Array{Any,1}
function GroupAlgebraElement(coordinates::Vector{T},
product_matrix::Array{Int,2})
length(coordinates) == size(product_matrix,1) ||
throw(ArgumentError("Matrix has to represent products of basis
elements"))
size(product_matrix, 1) == size(product_matrix, 2) ||
throw(ArgumentError("Product matrix has to be square"))
# length(coordinates) == length(basis) || throw(ArgumentError("Coordinates must be given in the given basis"))
# new(coordinates, product_matrix, basis)
new(coordinates, product_matrix)
end
end
# GroupAlgebraElement(c,pm,b) = GroupAlgebraElement(c,pm)
GroupAlgebraElement{T}(c::Vector{T},pm) = GroupAlgebraElement{T}(c,pm)
convert{T<:Number}(::Type{T}, X::GroupAlgebraElement) =
GroupAlgebraElement(convert(Vector{T}, X.coordinates), X.product_matrix)
show{T}(io::IO, X::GroupAlgebraElement{T}) = print(io,
"Element of Group Algebra over ", T, "\n", X.coordinates)
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.coordinates == Y.coordinates || 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(DomainError(
"Elements don't seem to belong to the same Group Algebra!"))
return GroupAlgebraElement(X.coordinates+Y.coordinates, 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.coordinates, Y.coordinates)...),
X.product_matrix)
end
(+)(X::GroupAlgebraElement, Y::GroupAlgebraElement) = add(X,Y)
(-)(X::GroupAlgebraElement) = GroupAlgebraElement(-X.coordinates, X.product_matrix)
(-)(X::GroupAlgebraElement, Y::GroupAlgebraElement) = add(X,-Y)
function group_star_multiplication{T<:Number}(X::GroupAlgebraElement{T},
Y::GroupAlgebraElement{T})
X.product_matrix == Y.product_matrix || DomainError(
"Elements don't seem to belong to the same Group Algebra!")
result = zeros(X.coordinates)
for (i,x) in enumerate(X.coordinates), (j,y) in enumerate(Y.coordinates)
index = X.product_matrix[i,j]
if index != 0
result[index]+= x*y
end
end
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.coordinates, X.product_matrix)
function scalar_multiplication{T<:Number, S<:Number}(a::T,
X::GroupAlgebraElement{S})
if T!=S
warn("Scalars and coefficients ring are not the same! Trying to promote...")
end
return GroupAlgebraElement(a*X.coordinates, 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.coordinates//a, X.product_matrix)
//{T<:Rational, S<:Integer}(X::GroupAlgebraElement{T}, a::S) =
X//convert(T,a)
length(X::GroupAlgebraElement) = length(X.coordinates)
size(X::GroupAlgebraElement) = size(X.coordinates)
norm(X::GroupAlgebraElement, p=2) = norm(X.coordinates, p)
end

117
Matrix_Constraints.g Normal file
View File

@ -0,0 +1,117 @@
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;
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;
USupport := function(x)
return Unique(Support(x));
end;
IsSupportedOn := function(basis, elt)
local elt_supp, x;
elt_supp := USupport(elt);
for x in elt_supp do
if not x in basis then
return x;
fi;
od;
return true;
end;
SDPGenerateAll := function(G, S, basis, name)
local QG, emb, delta, delta_sq, delta_vec, delta_sq_vec, product_constr;
QG := GroupRing(Rationals, G);;
emb := Embedding(G,QG);;
delta := Laplacian(G, S);;
delta_sq := delta^2;;
if not IsSupportedOn(basis, delta_sq) then
# Print("delta_sq is not supported on basis\n");
return fail;
else
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;;

70
property(T).jl Normal file
View File

@ -0,0 +1,70 @@
using JuMP
function read_GAP_raw_list(filename::String)
return eval(parse(String(read(filename))))
end
function create_product_matrix(matrix_constraints)
l = length(matrix_constraints)
product_matrix = zeros(Int, (l, l))
for (index, pairs) in enumerate(matrix_constraints)
for (i,j) in pairs
product_matrix[i,j] = index
end
end
return product_matrix
end
function create_sparse_product_matrix(matrix_constraints)
row_indices = Int[]
column_indices = Int[]
values = Int[]
for (index, pairs) in enumerate(matrix_constraints)
for (i,j) in pairs
push!(row_indices, i)
push!(column_indices, j)
push!(values, index)
end
end
sparse_product_matrix = sparse(row_indices, column_indices, values)
return sparse_product_matrix
end
function create_SDP_problem(matrix_constraints,
Δ²::GroupAlgebraElement, Δ::GroupAlgebraElement)
N = length(Δ)
@assert length(Δ) == length(Δ²)
@assert length(Δ) == length(matrix_constraints)
m = Model();
@variable(m, A[1:N, 1:N], SDP)
@SDconstraint(m, A >= zeros(size(A)))
@variable(m, κ >= 0.0)
@objective(m, Max, κ)
for (pairs, δ², δ) in zip(matrix_constraints, Δ².coordinates, Δ.coordinates)
@constraint(m, sum(A[i,j] for (i,j) in pairs) == δ² - κ*δ)
end
return m
end
function resulting_SOS{T<:Number, S<:Number}(sqrt_matrix::Array{T,2}, elt::GroupAlgebraElement{S})
zzz = zeros(T, size(sqrt_matrix)[1])
result::GroupAlgebraElement{T} = GroupAlgebraElement(zzz, elt.product_matrix)
for i in 1:length(result)
new_base = GroupAlgebraElement(sqrt_matrix[:,i], elt.product_matrix)
result += new_base*new_base
end
return result
end
function correct_to_augmentation_ideal{T<:Rational}(sqrt_matrix::Array{T,2})
sqrt_corrected = similar(sqrt_matrix)
l = size(sqrt_matrix,2)
for i in 1:l
col = view(sqrt_matrix,:,i)
sqrt_corrected[:,i] = col - sum(col)//l
# @assert sum(sqrt_corrected[:,i]) == 0
end
return sqrt_corrected
end