diff --git a/src/1712.07167.jl b/src/1712.07167.jl index 6e28349..36abe6c 100644 --- a/src/1712.07167.jl +++ b/src/1712.07167.jl @@ -89,7 +89,7 @@ function approximate_by_SOS(sett::Naive, isdir(fullpath(sett)) || mkpath(fullpath(sett)) @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 "Logging solver's progress into $solverlog" @@ -141,7 +141,7 @@ function approximate_by_SOS(sett::Symmetrized, orbit_data = decimate(orbit_data) @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 "Logging solver's progress into $solverlog" diff --git a/src/sos_sdps.jl b/src/sos_sdps.jl index 29ad03d..2853c8f 100644 --- a/src/sos_sdps.jl +++ b/src/sos_sdps.jl @@ -43,7 +43,7 @@ end ############################################################################### function SOS_problem_dual(elt::GroupRingElem, order_unit::GroupRingElem; - upper_bound::Float64=Inf) + lower_bound::Float64=Inf) @assert parent(elt) == parent(order_unit) 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, psd, [y[i] for i in RG.pm] in PSDCone()) - if !isinf(upper_bound) + if !isinf(lower_bound) @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 else @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]) # SP = Symmetric(P) - JuMP.@constraint(m, sdp, P in PSDCone()) + JuMP.@SDconstraint(m, sdp, P >= 0) if iszero(aug(X)) && iszero(aug(orderunit)) JuMP.@constraint(m, augmentation, sum(P) == 0) @@ -97,15 +98,13 @@ function SOS_problem_primal(X::GroupRingElem, orderunit::GroupRingElem; return m end -const SOS_problem = SOS_problem_primal - ############################################################################### # # 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) m = JuMP.Model(); @@ -113,7 +112,7 @@ function SOS_problem(X::GroupRingElem, orderunit::GroupRingElem, data::OrbitData for (k,s) in enumerate(Ns) Ps[k] = JuMP.@variable(m, [1:s, 1:s]) - JuMP.@constraint(m, Ps[k] in PSDCone()) + JuMP.@SDconstraint(m, Ps[k] >= 0) end if upper_bound < Inf diff --git a/test/1812.03456.jl b/test/1812.03456.jl index aca1905..78acc7e 100644 --- a/test/1812.03456.jl +++ b/test/1812.03456.jl @@ -100,7 +100,7 @@ end end 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); Base.Libc.flush_cstdio()