kalmar
/
GroupsWithPropertyT
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

last minute functional changes

This commit is contained in:
kalmarek 2018-11-24 13:58:49 +01:00
parent 9fc23c43d7
commit 8666cc6b2c
1 changed files with 39 additions and 33 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

38
positivity/check_positivity.jl
View File

@ -129,9 +129,13 @@ orbit_data = PropertyT.decimate(orbit_data);
elt = adj+3op; elt = adj+3op;
SDP_problem, varλ, varP = PropertyT.SOS_problem(elt, Δ, orbit_data; upper_bound=sett.upper_bound) const SOLUTION_FILE = PropertyT.filename(sett, :solution)
begin if !isfile(SOLUTION_FILE)
SDP_problem, varλ, varP = PropertyT.SOS_problem(elt, Δ, orbit_data; upper_bound=sett.upper_bound)
begin
using SCS using SCS
scs_solver = SCS.SCSSolver(linearsolver=SCS.Direct, scs_solver = SCS.SCSSolver(linearsolver=SCS.Direct,
eps=sett.tol, eps=sett.tol,
@ -140,15 +144,15 @@ begin
acceleration_lookback=1) acceleration_lookback=1)
JuMP.setsolver(SDP_problem, scs_solver) JuMP.setsolver(SDP_problem, scs_solver)
end end
λ = Ps = nothing λ = Ps = nothing
ws = PropertyT.warmstart(sett) ws = PropertyT.warmstart(sett)
# using ProgressMeter # using ProgressMeter
# @showprogress "Running SDP optimization step... " for i in 1:args["repetitions"] # @showprogress "Running SDP optimization step... " for i in 1:args["repetitions"]
while true while true
λ, Ps, ws = PropertyT.solve(PropertyT.filename(sett, :solverlog), λ, Ps, ws = PropertyT.solve(PropertyT.filename(sett, :solverlog),
SDP_problem, varλ, varP, ws); SDP_problem, varλ, varP, ws);
@ -158,11 +162,15 @@ while true
else else
warn("No valid solution was saved!") warn("No valid solution was saved!")
end end
end
info("Reconstructing P...")
@time P = PropertyT.reconstruct(Ps, orbit_data);
save(SOLUTION_FILE, "λ", λ, "P", P)
end end
info("Reconstructing P...") P, λ = load(SOLUTION_FILE, "P", "λ")
@time P = PropertyT.reconstruct(Ps, orbit_data); @show λ;
save(PropertyT.filename(sett, :solution), "λ", λ, "P", P)
@time const Q = real(sqrtm(P)); @time const Q = real(sqrtm(P));
@ -172,14 +180,12 @@ function SOS_residual(eoi::GroupRingElem, Q::Matrix)
return eoi - sos return eoi - sos
end end
addprocs(4) info("Floating Point arithmetic:")
@everywhere using PropertyT
@show λ;
EOI = elt - λ*Δ EOI = elt - λ*Δ
b = SOS_residual(EOI, Q) b = SOS_residual(EOI, Q)
@show norm(b, 1); @show norm(b, 1);
info("Interval arithmetic:")
using IntervalArithmetic using IntervalArithmetic
Qint = PropertyT.augIdproj(Q); Qint = PropertyT.augIdproj(Q);
@assert all([zero(eltype(Q)) in sum(view(Qint, :, i)) for i in 1:size(Qint, 2)]) @assert all([zero(eltype(Q)) in sum(view(Qint, :, i)) for i in 1:size(Qint, 2)])
@ -190,4 +196,4 @@ Q_int = PropertyT.augIdproj(Q);
b_int = SOS_residual(EOI_int, Q_int) b_int = SOS_residual(EOI_int, Q_int)
@show norm(b_int, 1); @show norm(b_int, 1);
@interval(λ) - 2^2*norm(b_int,1) info("λ is certified to be > ", (@interval(λ) - 2^2*norm(b_int,1)).lo)