diff --git a/MCG.jl b/MCG.jl index fba925c..74e9394 100644 --- a/MCG.jl +++ b/MCG.jl @@ -1,110 +1,187 @@ using ArgParse -using JLD - -using Nemo -import SCS.SCSSolver -using PropertyT - -using Groups - -function cpuinfo_physicalcores() - maxcore = -1 - for line in eachline("/proc/cpuinfo") - if startswith(line, "core id") - maxcore = max(maxcore, parse(Int, split(line, ':')[2])) - end - end - maxcore < 0 && error("failure to read core ids from /proc/cpuinfo") - return maxcore + 1 -end function parse_commandline() args = ArgParseSettings() - @add_arg_table args begin - "--tol" - help = "set numerical tolerance for the SDP solver" - arg_type = Float64 - default = 1e-14 - "--iterations" - help = "set maximal number of iterations for the SDP solver (default: 20000)" - arg_type = Int - default = 60000 - "--upper-bound" - help = "Set an upper bound for the spectral gap" - arg_type = Float64 - default = Inf - "--cpus" - help = "Set number of cpus used by solver (default: auto)" - arg_type = Int - required = false - "-N" - help = "Consider mapping class group of surface of genus N" - arg_type = Int - default = 2 - "--radius" - help = "Radius of ball B_r(e,S) to find solution over" - arg_type = Int - default = 4 - "--warmstart" - help = "Use warmstart.jl as the initial guess for SCS" - action = :store_true - end + @add_arg_table args begin + "--tol" + help = "set numerical tolerance for the SDP solver" + arg_type = Float64 + default = 1e-14 + "--iterations" + help = "set maximal number of iterations for the SDP solver (default: 20000)" + arg_type = Int + default = 60000 + "--upper-bound" + help = "Set an upper bound for the spectral gap" + arg_type = Float64 + default = Inf + "--cpus" + help = "Set number of cpus used by solver (default: auto)" + arg_type = Int + required = false + "-N" + help = "Consider mapping class group of surface of genus N" + arg_type = Int + default = 2 + "--radius" + help = "Radius of ball B_r(e,S) to find solution over" + arg_type = Int + default = 4 + "--warmstart" + help = "Use warmstart.jl as the initial guess for SCS" + action = :store_true + end - return parse_args(args) + return parse_args(args) end +const PARSEDARGS = parse_commandline() + +include("CPUselect.jl") +# set_parallel_mthread(PARSEDARGS, workers=true) include("FPGroups_GAP.jl") -include("CPUselect.jl") -function main() +module MCGrps - parsed_args = parse_commandline() - set_parallel_mthread(parsed_args) +using Groups +using Nemo - tol = parsed_args["tol"] - iterations = parsed_args["iterations"] - upper_bound = parsed_args["upper-bound"] - radius = parsed_args["radius"] - N = parsed_args["N"] +Comm(x,y) = x*y*x^-1*y^-1 - prefix = "MCG($N)" - name = "$(prefix)" +function Group(N::Int) + if N == 2 + MCGroup = Groups.FPGroup(["a1","a2","a3","a4","a5"]); + S = Nemo.gens(MCGroup) - isdir(name) || mkdir(name) + N = length(S) + A = prod(reverse(S))*prod(S) - prepare_pm_delta(prefix, name, radius) + relations = [ + [Comm(S[i], S[j]) for i in 1:N for j in 1:N if abs(i-j) > 1]..., + [S[i]*S[i+1]*S[i]*inv(S[i+1]*S[i]*S[i+1]) for i in 1:N-1]..., + (S[1]*S[2]*S[3])^4*inv(S[5])^2, + Comm(A, S[1]), + A^2 + ] - logger = PropertyT.setup_logging(name) + relations = [relations; [inv(rel) for rel in relations]] + Groups.add_rels!(MCGroup, Dict(rel => MCGroup() for rel in relations)) + return MCGroup - info(logger, "Group: $name") - info(logger, "Iterations: $iterations") - info(logger, "Precision: $tol") - info(logger, "Upper bound: $upper_bound") + elseif N < 2 + throw("Genus must be at least 2!") + end - MCGroup = Groups.FPGroup(["a1","a2","a3","a4", "a5"]); - S = Nemo.gens(MCGroup) - Comm(x,y) = x*y*x^-1*y^-1 - k = length(S) + MCGroup = Groups.FPGroup(["a$i" for i in 0:2N]) + S = Nemo.gens(MCGroup) - relations = [[Comm(S[i], S[j]) for i in 1:k for j in 1:k if abs(i-j) > 1]..., - [S[i]*S[i+1]*S[i]*inv(S[i+1]*S[i]*S[i+1]) for i in 1:k-1]..., - (S[1]*S[2]*S[3])^4*inv(S[5])^5, - Comm(prod(reverse(S))*prod(S), S[1]), - (prod(reverse(S))*prod(S))^2 - ]; + a0 = S[1] + A = S[2:end] + k = length(A) - relations = [relations; [inv(rel) for rel in relations]] - Groups.add_rels!(MCGroup, Dict(rel => MCGroup() for rel in relations)) + relations = [ + [Comm(A[i], A[j]) for i in 1:k for j in 1:k if abs(i-j) > 1]..., + [Comm(a0, A[i]) for i in 1:k if i != 4]..., + [A[i]*A[(i+1)]*A[i]*inv(A[i+1]*A[i]*A[i+1]) for i in 1:k-1]..., + A[4]*a0*A[4]*inv(a0*A[4]*a0) + ] - S = gens(MCGroup) - S = unique([S; [inv(s) for s in S]]) - Id = MCGroup() - - solver = SCSSolver(eps=tol, max_iters=iterations, linearsolver=SCS.Direct, alpha=1.9, acceleration_lookback=1) + # 3-chain relation + c = prod(reverse(A[1:4]))*prod(A[1:4]) + b0 = c*a0*inv(c) + push!(relations, (A[1]*A[2]*A[3])^4*inv(a0*b0)) - @time PropertyT.check_property_T(name, S, Id, solver, upper_bound, tol, radius) - return 0 + # Lantern relation + b1 = inv(A[4]*A[5]*A[3]*A[4])*a0*(A[4]*A[5]*A[3]*A[4]) + b2 = inv(A[2]*A[3]*A[1]*A[2])*b1*(A[2]*A[3]*A[1]*A[2]) + u = inv(A[6]*A[5])*b1*(A[6]*A[5]) + x = prod(reverse(A[2:6]))*u*prod(inv.(A[1:4])) + b3 = x*a0*inv(x) + push!(relations, a0*b2*b1*inv(A[1]*A[3]*A[5]*b3)) + + # Hyperelliptic relation + X = prod(reverse(A))*prod(A) + + function n(i::Int, b=b0) + if i == 1 + return A[1] + elseif i == 2 + return b + else + return w(i-2)*n(i-2)*w(i-2) + end + end + + function w(i::Int) + (A[2i+4]*A[2i+3]*A[2i+2]* n(i+1))*(A[2i+1]*A[2i] *A[2i+2]*A[2i+1])* + (A[2i+3]*A[2i+2]*A[2i+4]*A[2i+3])*( n(i+1)*A[2i+2]*A[2i+1]*A[2i] ) + end + + push!(relations, X*n(N)*inv(n(N)*X)) + + relations = [relations; [inv(rel) for rel in relations]] + Groups.add_rels!(MCGroup, Dict(rel => MCGroup() for rel in relations)) + @show MCGroup + return MCGroup end -main() +############################################################################### +# +# Misc +# +############################################################################### + +function groupname(parsed_args) + N = parsed_args["N"] + return groupname(N), N +end + +groupname(N::Int) = "MCG$(N)" + +end #of module MCGrps + + + +using SCS.SCSSolver +using PropertyT + +function main(GROUP, parsed_args) + + radius = parsed_args["radius"] + tol = parsed_args["tol"] + iterations = parsed_args["iterations"] + upper_bound = parsed_args["upper-bound"] + warm = parsed_args["warmstart"] + + name, N = GROUP.groupname(parsed_args) + isdir(name) || mkdir(name) + + G = GROUP.Group(N) + S = Nemo.gens(G) + relations = [k*inv(v) for (k,v) in G.rels] + prepare_pm_delta(name, GAP_groupcode(S, relations), radius) + + S = unique([S; [inv(s) for s in S]]) + + Id = G() + + logger = PropertyT.setup_logging(joinpath(name, "$(upper_bound)")) + + info(logger, "Group: $name") + info(logger, "Iterations: $iterations") + info(logger, "Precision: $tol") + info(logger, "Upper bound: $upper_bound") + info(logger, "Radius: $radius") + info(logger, G) + info(logger, "Symmetric generating set of size $(length(S))") + info(logger, "Threads: $(Threads.nthreads())") + info(logger, "Workers: $(workers())") + + solver = SCSSolver(eps=tol, max_iters=iterations, linearsolver=SCS.Direct, alpha=1.9, acceleration_lookback=1) + + PropertyT.check_property_T(name, S, Id, solver, upper_bound, tol, radius) + return 0 +end + +main(MCGrps, PARSEDARGS)