From 70c160411da5871c51e7b39dc955300da86ceeab Mon Sep 17 00:00:00 2001 From: kalmarek Date: Fri, 12 Apr 2019 23:18:48 +0200 Subject: [PATCH] cosmetic changes (mostly empty lines) --- src/1712.07167.jl | 66 +++++++++++++++++++++++--------------------- src/checksolution.jl | 4 +-- src/laplacians.jl | 2 +- src/orbitdata.jl | 6 ++-- src/sos_sdps.jl | 26 ++++++++--------- 5 files changed, 54 insertions(+), 50 deletions(-) diff --git a/src/1712.07167.jl b/src/1712.07167.jl index 093f8cd..dad4a61 100644 --- a/src/1712.07167.jl +++ b/src/1712.07167.jl @@ -14,7 +14,7 @@ struct Naive{El} <: Settings S::Vector{El} radius::Int upper_bound::Float64 - + solver::JuMP.OptimizerFactory warmstart::Bool end @@ -26,7 +26,7 @@ struct Symmetrized{El} <: Settings autS::Group radius::Int upper_bound::Float64 - + solver::JuMP.OptimizerFactory warmstart::Bool end @@ -87,39 +87,43 @@ function warmstart(sett::Settings) return ws end -function approximate_by_SOS(sett::Naive, +function approximate_by_SOS(sett::Naive, elt::GroupRingElem, orderunit::GroupRingElem; solverlog=tempname()*".log") - + + isdir(fullpath(sett)) || mkpath(fullpath(sett)) + @info "Creating SDP problem..." SDP_problem = SOS_problem(elt, orderunit, upper_bound=sett.upper_bound) @info Base.repr(SDP_problem) - + @info "Logging solver's progress into $solverlog" ws = warmstart(sett) @time status, ws = PropertyT.solve(solverlog, SDP_problem, sett.solver, ws) @info "Optimization finished:" status - + P = value.(SDP_problem[:P]) λ = value(SDP_problem[:λ]) - + if any(isnan.(P)) @warn "The solution seems to contain NaNs. Not overriding warmstart.jld" else save(filename(sett, :warmstart), "warmstart", (ws.primal, ws.dual, ws.slack), "P", P, "λ", λ) end - + save(filename(sett, :warmstart, date=true), "warmstart", (ws.primal, ws.dual, ws.slack), "P", P, "λ", λ) return λ, P end -function approximate_by_SOS(sett::Symmetrized, +function approximate_by_SOS(sett::Symmetrized, elt::GroupRingElem, orderunit::GroupRingElem; solverlog=tempname()*".log") - + + isdir(fullpath(sett)) || mkpath(fullpath(sett)) + if !isfile(filename(sett, :OrbitData)) isdefined(parent(orderunit), :basis) || throw("You need to define basis of Group Ring to compute orbit decomposition!") orbit_data = OrbitData(parent(orderunit), sett.autS) @@ -131,30 +135,30 @@ function approximate_by_SOS(sett::Symmetrized, @info "Creating SDP problem..." SDP_problem, varP = SOS_problem(elt, orderunit, orbit_data, upper_bound=sett.upper_bound) @info Base.repr(SDP_problem) - + @info "Logging solver's progress into $solverlog" - + ws = warmstart(sett) @time status, ws = PropertyT.solve(solverlog, SDP_problem, sett.solver, ws) @info "Optimization finished:" status - + λ = value(SDP_problem[:λ]) Ps = [value.(P) for P in varP] - + if any(any(isnan.(P)) for P in Ps) @warn "The solution seems to contain NaNs. Not overriding warmstart.jld" else save(filename(sett, :warmstart), "warmstart", (ws.primal, ws.dual, ws.slack), "Ps", Ps, "λ", λ) end - - save(filename(sett, :warmstart, date=true), + + save(filename(sett, :warmstart, date=true), "warmstart", (ws.primal, ws.dual, ws.slack), "Ps", Ps, "λ", λ) - + @info "Reconstructing P..." @time P = reconstruct(Ps, orbit_data) return λ, P -end +end ############################################################################### # @@ -167,21 +171,21 @@ function certify_SOS_decomposition(elt::GroupRingElem, orderunit::GroupRingElem, separator = "-"^76 @info "$separator\nChecking in floating-point arithmetic..." λ eoi = elt - λ*orderunit - + @info("Computing sum of squares decomposition...") @time residual = eoi - compute_SOS(parent(eoi), augIdproj(Q)) L1_norm = norm(residual,1) - certified_λ = λ - 2.0^(2ceil(log2(R)))*L1_norm - + floatingpoint_λ = λ - 2.0^(2ceil(log2(R)))*L1_norm + info_strs = ["Numerical metrics of the obtained SOS:", "ɛ(elt - λu - ∑ξᵢ*ξᵢ) ≈ $(aug(residual))", "‖elt - λu - ∑ξᵢ*ξᵢ‖₁ ≈ $(L1_norm)", "Floating point (NOT certified) λ ≈"] - @info join(info_strs, "\n") certified_λ + @info join(info_strs, "\n") floatingpoint_λ - if certified_λ ≤ 0 - return certified_λ + if floatingpoint_λ ≤ 0 + return floatingpoint_λ end λ = @interval(λ) @@ -190,18 +194,18 @@ function certify_SOS_decomposition(elt::GroupRingElem, orderunit::GroupRingElem, "λ ∈ $λ"] @info(join(info_strs, "\n")) eoi = elt - λ*orderunit - + @info("Projecting columns of Q to the augmentation ideal...") @time Q, check = augIdproj(Interval, Q) @info "Checking that sum of every column contains 0.0..." check_augmented=check check || @warn("The following numbers are meaningless!") - + @info("Computing sum of squares decomposition...") @time residual = eoi - compute_SOS(parent(eoi), Q) - + L1_norm = norm(residual,1) certified_λ = λ - 2.0^(2ceil(log2(R)))*L1_norm - + info_strs = ["Numerical metrics of the obtained SOS:", "ɛ(elt - λu - ∑ξᵢ*ξᵢ) ∈ $(aug(residual))", "‖elt - λu - ∑ξᵢ*ξᵢ‖₁ ∈ $(L1_norm)", @@ -262,7 +266,7 @@ end function spectral_gap(sett::Settings) fp = PropertyT.fullpath(sett) isdir(fp) || mkpath(fp) - + if isfile(filename(sett,:Δ)) # cached @info "Loading precomputed Δ..." @@ -283,7 +287,7 @@ function spectral_gap(sett::Settings) λ < 0 && @warn "Solver did not produce a valid solution!" end - + info_strs = ["Numerical metrics of matrix solution:", "sum(P) = $(sum(P))", "maximum(P) = $(maximum(P))", @@ -295,6 +299,6 @@ function spectral_gap(sett::Settings) @time Q = real(sqrt(Symmetric( (P.+ P')./2 ))) certified_sgap = spectral_gap(Δ, λ, Q, R=sett.radius) - + return certified_sgap end diff --git a/src/checksolution.jl b/src/checksolution.jl index c599482..e99fca2 100644 --- a/src/checksolution.jl +++ b/src/checksolution.jl @@ -39,11 +39,11 @@ end function compute_SOS_square(RG::GroupRing, Q::AbstractMatrix{<:Real}) result = zeros(eltype(Q), maximum(RG.pm)); - + for i in 1:size(Q,2) GroupRings.fmac!(result, view(Q,:,i), view(Q,:,i), RG.pm) end - + return GroupRingElem(result, RG) end diff --git a/src/laplacians.jl b/src/laplacians.jl index 252fcdc..d30bc13 100644 --- a/src/laplacians.jl +++ b/src/laplacians.jl @@ -59,7 +59,7 @@ end function loadGRElem(fname::String, G::Group) pm = load(fname, "pm") RG = GroupRing(G, pm) - return loadGRElem(fname, RG) + return loadGRElem(fname, RG) end function loadGRElem(fname::String) diff --git a/src/orbitdata.jl b/src/orbitdata.jl index 6cfe73f..4a8b61f 100644 --- a/src/orbitdata.jl +++ b/src/orbitdata.jl @@ -181,7 +181,7 @@ end function (p::perm)(A::GroupRingElem) RG = parent(A) result = zero(RG, eltype(A.coeffs)) - + for (idx, c) in enumerate(A.coeffs) if c!= zero(eltype(A.coeffs)) result[p(RG.basis[idx])] = c @@ -218,8 +218,8 @@ import Base.* """ function *(x::AbstractAlgebra.MatElem, P::Generic.perm) z = similar(x) - m = rows(x) - n = cols(x) + m = nrows(x) + n = ncols(x) for i = 1:m for j = 1:n z[i, j] = x[i,P[j]] diff --git a/src/sos_sdps.jl b/src/sos_sdps.jl index eb27bc7..5c81d6b 100644 --- a/src/sos_sdps.jl +++ b/src/sos_sdps.jl @@ -55,14 +55,14 @@ function SOS_problem(X::GroupRingElem, orderunit::GroupRingElem; upper_bound::Fl else λ = JuMP.@variable(m, λ) end - + cnstrs = constraints(parent(X).pm) for (constraint_indices, x, u) in zip(cnstrs, X.coeffs, orderunit.coeffs) JuMP.@constraint(m, x - λ*u == sum(P[constraint_indices])) end - + JuMP.@objective(m, Max, λ) - + return m end @@ -88,12 +88,12 @@ function SOS_problem(X::GroupRingElem, orderunit::GroupRingElem, data::OrbitData else λ = JuMP.@variable(m, λ) end - + @info "Adding $(length(data.orbits)) constraints..." @time addconstraints!(m, Ps, X, orderunit, data) JuMP.@objective(m, Max, λ) - + return m, Ps end @@ -121,9 +121,9 @@ function addconstraints!(m::JuMP.Model, for (t, orbit) in enumerate(data.orbits) orbit_constraint!(orb_cnstr, cnstrs, orbit) constraintLHS!(M, orb_cnstr, data.Uπs, UπsT, data.dims) - + x, u = X_orb[t], orderunit_orb[t] - + JuMP.@constraints m begin x - λ*u == sum(dot(M[π], P[π]) for π in eachindex(data.Uπs)) end @@ -142,7 +142,7 @@ function reconstruct(Ps::Vector{M}, lU = length(Uπs) transfP = [dims[π].*Uπs[π]*Ps[π]*Uπs[π]' for π in 1:lU] tmp = [zeros(Float64, size(first(transfP))) for _ in 1:lU] - + Threads.@threads for π in 1:lU tmp[π] = perm_avg(tmp[π], transfP[π], values(preps)) end @@ -172,7 +172,7 @@ end # ############################################################################### -function warmstart_scs!(m::JuMP.Model, warmstart) +function setwarmstart_scs!(m::JuMP.Model, warmstart) solver_name(m) == "SCS" || throw("warmstarting defined only for SCS!") primal, dual, slack = warmstart m.moi_backend.optimizer.model.optimizer.sol.primal = primal @@ -190,14 +190,14 @@ function getwarmstart_scs(m::JuMP.Model) ) return warmstart end - + function solve(m::JuMP.Model, with_optimizer::JuMP.OptimizerFactory, warmstart=nothing) - + set_optimizer(m, with_optimizer) MOIU.attach_optimizer(m) - + if warmstart != nothing - warmstart_scs!(m, warmstart) + setwarmstart_scs!(m, warmstart) end optimize!(m)