1
0
mirror of https://github.com/kalmarek/PropertyT.jl.git synced 2024-11-26 00:55:27 +01:00

use @SDconstraint and SOS_problem_(primal|dual)

This commit is contained in:
kalmarek 2020-10-17 02:11:34 +02:00
parent 7ab971eb76
commit eeecc35232
No known key found for this signature in database
GPG Key ID: 8BF1A3855328FC15
3 changed files with 10 additions and 11 deletions

View File

@ -89,7 +89,7 @@ function approximate_by_SOS(sett::Naive,
isdir(fullpath(sett)) || mkpath(fullpath(sett)) isdir(fullpath(sett)) || mkpath(fullpath(sett))
@info "Creating SDP problem..." @info "Creating SDP problem..."
SDP_problem = SOS_problem(elt, orderunit, upper_bound=sett.upper_bound) SDP_problem = SOS_problem_primal(elt, orderunit, upper_bound=sett.upper_bound)
@info Base.repr(SDP_problem) @info Base.repr(SDP_problem)
@info "Logging solver's progress into $solverlog" @info "Logging solver's progress into $solverlog"
@ -141,7 +141,7 @@ function approximate_by_SOS(sett::Symmetrized,
orbit_data = decimate(orbit_data) orbit_data = decimate(orbit_data)
@info "Creating SDP problem..." @info "Creating SDP problem..."
SDP_problem, varP = SOS_problem(elt, orderunit, orbit_data, upper_bound=sett.upper_bound) SDP_problem, varP = SOS_problem_primal(elt, orderunit, orbit_data, upper_bound=sett.upper_bound)
@info Base.repr(SDP_problem) @info Base.repr(SDP_problem)
@info "Logging solver's progress into $solverlog" @info "Logging solver's progress into $solverlog"

View File

@ -43,7 +43,7 @@ end
############################################################################### ###############################################################################
function SOS_problem_dual(elt::GroupRingElem, order_unit::GroupRingElem; function SOS_problem_dual(elt::GroupRingElem, order_unit::GroupRingElem;
upper_bound::Float64=Inf) lower_bound::Float64=Inf)
@assert parent(elt) == parent(order_unit) @assert parent(elt) == parent(order_unit)
RG = parent(elt) RG = parent(elt)
@ -54,9 +54,10 @@ function SOS_problem_dual(elt::GroupRingElem, order_unit::GroupRingElem;
@constraint(m, λ_dual, dot(order_unit.coeffs, y) == 1) @constraint(m, λ_dual, dot(order_unit.coeffs, y) == 1)
@constraint(m, psd, [y[i] for i in RG.pm] in PSDCone()) @constraint(m, psd, [y[i] for i in RG.pm] in PSDCone())
if !isinf(upper_bound) if !isinf(lower_bound)
@variable(m, λ_ub_dual) @variable(m, λ_ub_dual)
expr = dot(elt.coeffs, y) + upper_bound*λ_ub_dual expr = dot(elt.coeffs, y) + lower_bound*λ_ub_dual
# @constraint m expr >= lower_bound
@objective m Min expr @objective m Min expr
else else
@objective m Min dot(elt.coeffs, y) @objective m Min dot(elt.coeffs, y)
@ -73,7 +74,7 @@ function SOS_problem_primal(X::GroupRingElem, orderunit::GroupRingElem;
JuMP.@variable(m, P[1:N, 1:N]) JuMP.@variable(m, P[1:N, 1:N])
# SP = Symmetric(P) # SP = Symmetric(P)
JuMP.@constraint(m, sdp, P in PSDCone()) JuMP.@SDconstraint(m, sdp, P >= 0)
if iszero(aug(X)) && iszero(aug(orderunit)) if iszero(aug(X)) && iszero(aug(orderunit))
JuMP.@constraint(m, augmentation, sum(P) == 0) JuMP.@constraint(m, augmentation, sum(P) == 0)
@ -97,15 +98,13 @@ function SOS_problem_primal(X::GroupRingElem, orderunit::GroupRingElem;
return m return m
end end
const SOS_problem = SOS_problem_primal
############################################################################### ###############################################################################
# #
# Symmetrized SDP # Symmetrized SDP
# #
############################################################################### ###############################################################################
function SOS_problem(X::GroupRingElem, orderunit::GroupRingElem, data::OrbitData; upper_bound::Float64=Inf) function SOS_problem_primal(X::GroupRingElem, orderunit::GroupRingElem, data::OrbitData; upper_bound::Float64=Inf)
Ns = size.(data.Uπs, 2) Ns = size.(data.Uπs, 2)
m = JuMP.Model(); m = JuMP.Model();
@ -113,7 +112,7 @@ function SOS_problem(X::GroupRingElem, orderunit::GroupRingElem, data::OrbitData
for (k,s) in enumerate(Ns) for (k,s) in enumerate(Ns)
Ps[k] = JuMP.@variable(m, [1:s, 1:s]) Ps[k] = JuMP.@variable(m, [1:s, 1:s])
JuMP.@constraint(m, Ps[k] in PSDCone()) JuMP.@SDconstraint(m, Ps[k] >= 0)
end end
if upper_bound < Inf if upper_bound < Inf

View File

@ -100,7 +100,7 @@ end
end end
function check_positivity(elt, Δ, orbit_data, upper_bound, warm=nothing; with_solver=with_SCS(20_000, accel=10)) function check_positivity(elt, Δ, orbit_data, upper_bound, warm=nothing; with_solver=with_SCS(20_000, accel=10))
SDP_problem, varP = PropertyT.SOS_problem(elt, Δ, orbit_data; upper_bound=upper_bound) SDP_problem, varP = PropertyT.SOS_problem_primal(elt, Δ, orbit_data; upper_bound=upper_bound)
status, warm = PropertyT.solve(SDP_problem, with_solver, warm); status, warm = PropertyT.solve(SDP_problem, with_solver, warm);
Base.Libc.flush_cstdio() Base.Libc.flush_cstdio()