From 50b7b574be1d3ea1717975b37d54a762246dcda6 Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Fri, 7 Sep 2018 19:49:29 -0400 Subject: [PATCH 01/10] WIP hook up to MOI; working on linear tests --- src/Alfonso.jl | 7 +- src/cone.jl | 19 +++ src/mathoptinterface.jl | 257 ++++++++++++++++++++++++++++------------ src/nativeinterface.jl | 2 +- test/moiexamples.jl | 43 ------- test/runtests.jl | 86 +++++++++----- 6 files changed, 261 insertions(+), 153 deletions(-) diff --git a/src/Alfonso.jl b/src/Alfonso.jl index 38c905639..54ea58bc5 100644 --- a/src/Alfonso.jl +++ b/src/Alfonso.jl @@ -4,11 +4,14 @@ module Alfonso using SparseArrays using LinearAlgebra - include("interpolation.jl") include("cone.jl") + include("interpolation.jl") for primcone in ["nonnegative", "sumofsquares", "secondorder", "exponential", "power", "rotatedsecondorder", "positivesemidefinite"] include(joinpath(@__DIR__, "primitivecones", primcone * ".jl")) end include("nativeinterface.jl") - # include("mathoptinterface.jl") + + import MathOptInterface + const MOI = MathOptInterface + include("mathoptinterface.jl") end diff --git a/src/cone.jl b/src/cone.jl index 3051aeb36..88c13bf85 100644 --- a/src/cone.jl +++ b/src/cone.jl @@ -9,6 +9,25 @@ abstract type PrimitiveCone end mutable struct Cone prms::Vector{PrimitiveCone} idxs::Vector{AbstractVector{Int}} + + function Cone(prms::Vector{PrimitiveCone}, idxs::Vector{AbstractVector{Int}}) + @assert length(prms) == length(idxs) + for k in eachindex(prms) + @assert dimension(prms[k]) == length(idxs[k]) + end + cone = new() + cone.prms = prms + cone.idxs = idxs + return cone + end +end +Cone() = Cone(PrimitiveCone[], AbstractVector{Int}[]) + +function addprimitivecone!(cone, prm, idx) + @assert dimension(prm) == length(idx) + push!(cone.prms, prm) + push!(cone.idxs, idx) + return cone end # calculate complexity parameter of the barrier (sum of the primitive cone barrier parameters) diff --git a/src/mathoptinterface.jl b/src/mathoptinterface.jl index b4779d6a8..7d6819f3b 100644 --- a/src/mathoptinterface.jl +++ b/src/mathoptinterface.jl @@ -1,9 +1,4 @@ -import MathOptInterface -const MOI = MathOptInterface -const MOIU = MOI.Utilities - - mutable struct Optimizer <: MOI.AbstractOptimizer alf::AlfonsoOpt end @@ -12,7 +7,7 @@ Optimizer() = Optimizer(AlfonsoOpt()) MOI.get(::Optimizer, ::MOI.SolverName) = "Alfonso" -MOI.isempty(opt::Optimizer) = (get_status(opt.alf) == :NotLoaded) +MOI.is_empty(opt::Optimizer) = (get_status(opt.alf) == :NotLoaded) MOI.empty!(opt::Optimizer) = Optimizer() # TODO empty the data and results, or just create a new one? keep options? MOI.supports(::Optimizer, ::MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}) = true @@ -20,27 +15,41 @@ MOI.supports(::Optimizer, ::MOI.ObjectiveSense) = true # TODO don't restrict to Float64 type const SupportedFuns = Union{ - # MOI.SingleVariable, - # MOI.ScalarAffineFunction{Float64}, + MOI.SingleVariable, + MOI.ScalarAffineFunction{Float64}, MOI.VectorOfVariables, MOI.VectorAffineFunction{Float64}, } + const SupportedSets = Union{ - # MOI.EqualTo{Float64}, - # MOI.GreaterThan{Float64}, - # MOI.LessThan{Float64}, + MOI.EqualTo{Float64}, + MOI.GreaterThan{Float64}, + MOI.LessThan{Float64}, MOI.Zeros, MOI.Nonnegatives, - # MOI.Nonpositives, + MOI.Nonpositives, + MOI.SecondOrderCone, + MOI.RotatedSecondOrderCone, + MOI.PositiveSemidefiniteConeTriangle, + MOI.ExponentialCone, + MOI.PowerCone, } -MOI.supportsconstraint(::Optimizer, ::Type{<:SupportedFuns}, ::Type{<:SupportedSets}) = true +MOI.supports_constraint(::Optimizer, ::Type{<:SupportedFuns}, ::Type{<:SupportedSets}) = true + +conefrommoi(s::MOI.Nonnegatives) = NonnegativeCone(MOI.dimension(s)) +conefrommoi(s::MOI.Nonpositives) = error("Nonpositive cone not handled yet") +conefrommoi(s::MOI.SecondOrderCone) = SecondOrderCone(MOI.dimension(s)) +conefrommoi(s::MOI.RotatedSecondOrderCone) = RotatedSecondOrderCone(MOI.dimension(s)) +conefrommoi(s::MOI.PositiveSemidefiniteConeTriangle) = PositiveSemidefiniteCone(MOI.dimension(s)) +conefrommoi(s::MOI.ExponentialCone) = ExponentialCone() +conefrommoi(s::MOI.PowerCone) = PowerCone(s.exponent) # build representation as min c'x s.t. Ax = b, x in K # TODO what if some variables x are in multiple cone constraints? -function MOI.copy!(opt::Optimizer, src::MOI.ModelLike; copynames=false, warnattributes=true) - @assert !copynames - idxmap = Dict{MOI.Index,Int}() +function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_attributes=true) + @assert !copy_names + idxmap = Dict{MOI.Index,MOI.Index}() # model # mattr_src = MOI.get(src, MOI.ListOfModelAttributesSet()) @@ -51,14 +60,14 @@ function MOI.copy!(opt::Optimizer, src::MOI.ModelLike; copynames=false, warnattr j = 0 for vj in MOI.get(src, MOI.ListOfVariableIndices()) # MOI.VariableIndex j += 1 - idxmap[vj] = j + idxmap[vj] = MOI.VariableIndex(j) end @assert j == n # objective function (Jc, Vc) = (Int[], Float64[]) - for t in (MOI.get(src, MOI.ObjectiveFunction())).terms - push!(Jc, idxmap[t.variable_index]) + for t in (MOI.get(src, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}())).terms + push!(Jc, idxmap[t.variable_index].value) push!(Vc, t.coefficient) end ismax = (MOI.get(src, MOI.ObjectiveSense()) == MOI.MaxSense) ? true : false @@ -66,79 +75,171 @@ function MOI.copy!(opt::Optimizer, src::MOI.ModelLike; copynames=false, warnattr Vc .*= -1 end - # constraints + getsrccons(F, S) = MOI.get(src, MOI.ListOfConstraintIndices{F,S}()) + getconfun(conidx) = MOI.get(src, MOI.ConstraintFunction(), conidx) + getconset(conidx) = MOI.get(src, MOI.ConstraintSet(), conidx) + + # pass over vector constraints to construct conic model data # TODO don't enforce Float64 type (IA, JA, VA) = (Int[], Int[], Float64[]) (Ib, Vb) = (Int[], Float64[]) - cones = ConeData[] - coneidxs = AbstractUnitRange[] + cone = Cone() i = 0 # MOI constraint objects m = 0 # rows of A - for (F, S) in MOI.get(src, MOI.ListOfConstraints()) - # fsattr_src = MOI.get(src, MOI.ListOfConstraintAttributesSet{F,S}) - for ci in MOI.get(src, MOI.ListOfConstraintIndices{F,S}()) # MOI.ConstraintIndex{F,S} - # TODO need to access row indices of constraints for get dual on equalities - i += 1 - idxmap[ci] = i - - fi = MOI.get(src, MOI.ConstraintFunction(), ci) - si = MOI.get(src, MOI.ConstraintSet(), ci) - - if isa(fi, MOI.AbstractScalarFunction) && isa(si, MOI.AbstractScalarSet) - error("scalar F and S not yet implemented") - elseif isa(fi, MOI.AbstractVectorFunction) && isa(si, MOI.AbstractVectorSet) - dim = MOI.dimension(si) - @assert MOI.output_dimension(fi) == dim - - if isa(fi, MOI.VectorOfVariables) - if isa(si, MOI.Zeros) - # variables equal to zero: error - error("variables cannot be set equal to zero") - else - # variables in cone: don't modify A and b - push!(cones, buildcone(si)) # TODO implement in barriers.jl - push!(coneidxs, [idxmap[vj] for vj in fi.variables]) - end - elseif isa(fi, MOI.VectorAffineFunction) - # vector affine function: modify A and b - append!(Ib, collect(m+1:m+dim)) - append!(Vb, -fi.constants) - - for vt in fi.terms - push!(IA, m+vt.output_index) - push!(JA, idxmap[vt.scalar_term.variable_index]) - push!(VA, vt.scalar_term.coefficient) - end - - if !isa(si, MOI.Zeros) - # add slack variables in new cone - append!(IA, collect(m+1:m+dim)) - append!(JA, collect(n+1:n+dim)) - append!(VA, fill(-1.0, dim)) - - push!(cones, buildcone(si)) # TODO implement in barriers.jl - push!(coneidxs, n+1:n+dim) - - n += dim - end - - m += dim - else - error("constraint with function $fi in set $si is not supported by Alfonso") - end - else - error("constraint with function $fi in set $si is not supported by Alfonso") - end + + # equality constraints + + for ci in getsrccons(MOI.SingleVariable, MOI.EqualTo{Float64}) + # TODO can preprocess out maybe + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, MOI.EqualTo{Float64}}(i) + + m += 1 + push!(Ib, m) + push!(Vb, getconset(ci).value) + push!(IA, m) + push!(JA, idxmap[getconfun(ci).variable].value) + push!(VA, 1.0) + end + + for ci in getsrccons(MOI.SingleVariable, MOI.Zeros) + # TODO can preprocess out maybe + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, MOI.Zeros}(i) + + m += 1 + push!(IA, m) + push!(JA, idxmap[getconfun(ci).variable].value) + push!(VA, 1.0) + end + + for ci in getsrccons(MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}}(i) + + fi = getconfun(ci) + m += 1 + push!(Ib, m) + push!(Vb, getconset(ci).value - fi.constant) + for vt in fi.terms + push!(IA, m) + push!(JA, idxmap[vt.variable_index].value) + push!(VA, vt.coefficient) + end + end + + for ci in getsrccons(MOI.ScalarAffineFunction{Float64}, MOI.Zeros) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, MOI.Zeros}(i) + + fi = getconfun(ci) + m += 1 + push!(Ib, m) + push!(Vb, -fi.constant) + for vt in fi.terms + push!(IA, m) + push!(JA, idxmap[vt.variable_index].value) + push!(VA, vt.coefficient) + end + end + + for ci in getsrccons(MOI.VectorOfVariables, MOI.Zeros) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, MOI.Zeros}(i) + + for vi in getconfun(ci).variables + m += 1 + push!(IA, m) + push!(JA, idxmap[vi].value) + push!(VA, 1.0) + end + end + + for ci in getsrccons(MOI.VectorAffineFunction{Float64}, MOI.Zeros) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, MOI.Zeros}(i) + + fi = getconfun(ci) + dim = MOI.output_dimension(fi) + append!(Ib, collect(m+1:m+dim)) + append!(Vb, -fi.constants) + for vt in fi.terms + push!(IA, m + vt.output_index) + push!(JA, idxmap[vt.scalar_term.variable_index].value) + push!(VA, vt.scalar_term.coefficient) + end + m += dim + end + + # linear inequality constraints + + nonnegvars = Int[] # for building up a single nonnegative cone + + for S in (MOI.GreaterThan{Float64}, MOI.LessThan{Float64}), ci in getsrccons(MOI.SingleVariable, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, S}(i) + + bval = ((S == MOI.GreaterThan{Float64}) ? getconset(ci).lower : getconset(ci).upper) + if iszero(bval) + # no need to change A,b + push!(nonnegvars, idxmap[getconfun(ci).variable].value) + else + # add auxiliary variable and equality constraint + m += 1 + push!(Ib, m) + push!(Vb, bval) + push!(IA, m) + push!(JA, idxmap[getconfun(ci).variable].value) + push!(VA, 1.0) + n += 1 + push!(IA, m) + push!(JA, n) + push!(VA, ((S == MOI.GreaterThan{Float64}) ? -1.0 : 1.0)) + push!(nonnegvars, n) end end + for S in (MOI.GreaterThan{Float64}, MOI.LessThan{Float64}), ci in getsrccons(MOI.ScalarAffineFunction{Float64}, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, S}(i) + + # add auxiliary variable and equality constraint + fi = getconfun(ci) + m += 1 + push!(Ib, m) + push!(Vb, ((S == MOI.GreaterThan{Float64}) ? getconset(ci).lower : getconset(ci).upper) - fi.constant) + for vt in fi.terms + push!(IA, m) + push!(JA, idxmap[vt.variable_index].value) + push!(VA, vt.coefficient) + end + n += 1 + push!(IA, m) + push!(JA, n) + push!(VA, ((S == MOI.GreaterThan{Float64}) ? -1.0 : 1.0)) + push!(nonnegvars, n) + end + + # MOI.VectorOfVariables, + # MOI.VectorAffineFunction{Float64}, + # MOI.Nonnegatives, + # MOI.Nonpositives, + + + # add single nonnegative cone + addprimitivecone!(cone, NonnegativeCone(length(nonnegvars)), nonnegvars) + + + # TODO iterate over types of S, from easiest to hardest, appending cones in order of toughness + # then SOC, RSOC, exp, power, sumofsquares, psd + A = dropzeros!(sparse(IA, JA, VA, m, n)) b = Vector(dropzeros!(sparsevec(Ib, Vb, m))) c = Vector(dropzeros!(sparsevec(Jc, Vc, n))) # load problem data through native interface - load_data!(opt.alf, A, b, c, cones, coneidxs) + load_data!(opt.alf, A, b, c, cone) return idxmap end diff --git a/src/nativeinterface.jl b/src/nativeinterface.jl index 45e09a319..e823bbcb7 100644 --- a/src/nativeinterface.jl +++ b/src/nativeinterface.jl @@ -84,7 +84,7 @@ mutable struct AlfonsoOpt end function AlfonsoOpt(; - verbose = false, + verbose = true, optimtol = 1e-6, maxiter = 1e3, predlinesearch = true, diff --git a/test/moiexamples.jl b/test/moiexamples.jl index b6e02588a..e69de29bb 100644 --- a/test/moiexamples.jl +++ b/test/moiexamples.jl @@ -1,43 +0,0 @@ - -# TODO -# @testset "envelope example" begin -# opt = build_envelope(2, 5, 1, 5, use_data=true, native=false) -# MOI.optimize!(opt) -# @test MOI.get(opt, MOI.TerminationStatus()) == MOI.Success -# @test MOI.get(opt, MOI.ObjectiveValue()) ≈ -25.502777 atol=1e-4 -# @test MOI.get(opt, MOI.ObjectiveBound()) ≈ -25.502777 atol=1e-4 -# end -# -# @testset "lp example" begin -# opt = build_lp(500, 1000, use_data=true, native=false) -# MOI.optimize!(opt) -# @test MOI.get(opt, MOI.TerminationStatus()) == MOI.Success -# @test MOI.get(opt, MOI.ObjectiveValue()) ≈ 2055.807 atol=1e-4 -# @test MOI.get(opt, MOI.ObjectiveBound()) ≈ 2055.807 atol=1e-4 -# end -# -# @testset "namedpoly examples" begin -# @testset "Goldstein-Price" begin -# opt = build_namedpoly(:goldsteinprice, 7, native=false) -# MOI.optimize!(opt) -# @test MOI.get(opt, MOI.TerminationStatus()) == MOI.Success -# @test MOI.get(opt, MOI.ObjectiveValue()) ≈ 3 atol=1e-4 -# @test MOI.get(opt, MOI.ObjectiveBound()) ≈ 3 atol=1e-4 -# end -# -# @testset "Robinson" begin -# opt = build_namedpoly(:robinson, 8, native=false) -# MOI.optimize!(opt) -# @test MOI.get(opt, MOI.TerminationStatus()) == MOI.Success -# @test MOI.get(opt, MOI.ObjectiveValue()) ≈ 0.814814 atol=1e-4 -# @test MOI.get(opt, MOI.ObjectiveBound()) ≈ 0.814814 atol=1e-4 -# end -# -# @testset "Lotka-Volterra" begin -# opt = build_namedpoly(:lotkavolterra, 3, native=false) -# MOI.optimize!(opt) -# @test MOI.get(opt, MOI.TerminationStatus()) == MOI.Success -# @test MOI.get(opt, MOI.ObjectiveValue()) ≈ -20.8 atol=1e-4 -# @test MOI.get(opt, MOI.ObjectiveBound()) ≈ -20.8 atol=1e-4 -# end -# end diff --git a/test/runtests.jl b/test/runtests.jl index e18934a2c..2a6f44d95 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2,33 +2,61 @@ using Alfonso using Test -verbflag = false # Alfonso verbose option - -# TODO interpolation tests - -# load optimizer builder functions from examples folder -egs_dir = joinpath(@__DIR__, "../examples") -include(joinpath(egs_dir, "envelope/envelope.jl")) -include(joinpath(egs_dir, "lp/lp.jl")) -include(joinpath(egs_dir, "namedpoly/namedpoly.jl")) - -# run native and MOI interfaces on examples -include(joinpath(@__DIR__, "nativeexamples.jl")) - -# using MathOptInterface -# const MOI = MathOptInterface -# const MOIT = MOI.Test -# const MOIB = MOI.Bridges -# const MOIU = MOI.Utilities - -# include(joinpath(@__DIR__, "moiexamples.jl")) - -# run MOI continuous linear and conic tests -# MOIU.@model AlfonsoModelData () () (Zeros, Nonnegatives) () () () (VectorOfVariables,) (VectorAffineFunction,) -# const optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Optimizer()) -# const config = MOIT.TestConfig(atol=1e-4, rtol=1e-4) - -# TODO use bridges e.g. GreaterThan constraints into Nonnegatives constraints -# @testset "Continuous linear problems" begin -# MOIT.contlineartest(MOIB.SplitInterval{Float64}(optimizer), config) +# verbflag = true # Alfonso verbose option + +# # TODO interpolation tests +# +# # load optimizer builder functions from examples folder +# egs_dir = joinpath(@__DIR__, "../examples") +# include(joinpath(egs_dir, "envelope/envelope.jl")) +# include(joinpath(egs_dir, "lp/lp.jl")) +# include(joinpath(egs_dir, "namedpoly/namedpoly.jl")) +# +# # run native and MOI interfaces on examples +# include(joinpath(@__DIR__, "nativeexamples.jl")) + + +import MathOptInterface +const MOI = MathOptInterface +const MOIT = MOI.Test +const MOIB = MOI.Bridges +const MOIU = MOI.Utilities + +MOIU.@model(AlfonsoModelData, + (), + (MOI.EqualTo, MOI.GreaterThan, MOI.LessThan), + (MOI.Zeros, MOI.Nonnegatives, MOI.Nonpositives, MOI.SecondOrderCone, MOI.RotatedSecondOrderCone, MOI.PositiveSemidefiniteConeTriangle, MOI.ExponentialCone, MOI.PowerCone), + (), + (MOI.SingleVariable,), + (MOI.ScalarAffineFunction,), + (MOI.VectorOfVariables,), + (MOI.VectorAffineFunction,), + ) + +const optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Optimizer()) + +const config = MOIT.TestConfig( + atol=1e-4, + rtol=1e-4, + solve=true, + query=false, + modify_lhs=false, + duals=false, + infeas_certificates=false, + ) + +@testset "Continuous linear problems" begin + MOIT.contlineartest(MOIB.SplitInterval{Float64}(optimizer), config) +end + +# @testset "Continuous conic problems" begin +# exclude = ["rootdet", "logdet"] +# MOIT.contconictest( +# MOIB.SquarePSD{Float64}( +# MOIB.GeoMean{Float64}( +# MOIB.LogDet{Float64}( +# MOIB.RootDet{Float64}( +# optimizer +# )))), +# config, exclude) # end From 0da5440e88706d07542db8837b323eb1f8963ab7 Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Sat, 8 Sep 2018 00:35:05 -0400 Subject: [PATCH 02/10] finish adding linear conic transformations --- src/cone.jl | 2 +- src/mathoptinterface.jl | 41 +++++++++++++++++++++++++++++++++++------ 2 files changed, 36 insertions(+), 7 deletions(-) diff --git a/src/cone.jl b/src/cone.jl index 88c13bf85..18158786e 100644 --- a/src/cone.jl +++ b/src/cone.jl @@ -23,7 +23,7 @@ mutable struct Cone end Cone() = Cone(PrimitiveCone[], AbstractVector{Int}[]) -function addprimitivecone!(cone, prm, idx) +function addprimitivecone!(cone::Cone, prm::PrimitiveCone, idx::AbstractVector{Int}) @assert dimension(prm) == length(idx) push!(cone.prms, prm) push!(cone.idxs, idx) diff --git a/src/mathoptinterface.jl b/src/mathoptinterface.jl index 7d6819f3b..e1ccf6d2c 100644 --- a/src/mathoptinterface.jl +++ b/src/mathoptinterface.jl @@ -181,16 +181,17 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, S}(i) bval = ((S == MOI.GreaterThan{Float64}) ? getconset(ci).lower : getconset(ci).upper) + vn = getconfun(ci).variable if iszero(bval) # no need to change A,b - push!(nonnegvars, idxmap[getconfun(ci).variable].value) + push!(nonnegvars, idxmap[vn].value) else # add auxiliary variable and equality constraint m += 1 push!(Ib, m) push!(Vb, bval) push!(IA, m) - push!(JA, idxmap[getconfun(ci).variable].value) + push!(JA, idxmap[vn].value) push!(VA, 1.0) n += 1 push!(IA, m) @@ -221,11 +222,39 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ push!(nonnegvars, n) end - # MOI.VectorOfVariables, - # MOI.VectorAffineFunction{Float64}, - # MOI.Nonnegatives, - # MOI.Nonpositives, + for S in (MOI.Nonnegatives, MOI.Nonpositives), ci in getsrccons(MOI.VectorOfVariables, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, S}(i) + + for vj in getconfun(ci).variables + col = idxmap[vj].value + m += 1 + push!(IA, m) + push!(JA, col) + push!(VA, ((S == MOI.Nonnegatives) ? 1.0 : -1.0)) + push!(nonnegvars, col) + end + end + + for S in (MOI.Nonnegatives, MOI.Nonpositives), ci in getsrccons(MOI.VectorAffineFunction{Float64}, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, S}(i) + fi = getconfun(ci) + dim = MOI.output_dimension(fi) + for vt in fi.terms + push!(IA, m + vt.output_index) + push!(JA, idxmap[vt.scalar_term.variable_index].value) + push!(VA, vt.scalar_term.coefficient) + end + append!(Ib, collect(m+1:m+dim)) + append!(Vb, -fi.constants) + append!(IA, collect(m+1:m+dim)) + append!(JA, collect(n+1:n+dim)) + append!(VA, fill(((S == MOI.Nonnegatives) ? 1.0 : -1.0), dim)) + m += dim + n += dim + end # add single nonnegative cone addprimitivecone!(cone, NonnegativeCone(length(nonnegvars)), nonnegvars) From d5bc54a02202f367cc0bdb7311c481ca35ad1c21 Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Sat, 8 Sep 2018 17:05:21 -0400 Subject: [PATCH 03/10] progress on getting linear tests to pass --- src/mathoptinterface.jl | 158 +++++++++++++++++- src/nativeinterface.jl | 15 +- test/runtests.jl | 362 +++++++++++++++++++++++++++++++++++++++- 3 files changed, 519 insertions(+), 16 deletions(-) diff --git a/src/mathoptinterface.jl b/src/mathoptinterface.jl index e1ccf6d2c..3c45f3e83 100644 --- a/src/mathoptinterface.jl +++ b/src/mathoptinterface.jl @@ -1,6 +1,13 @@ mutable struct Optimizer <: MOI.AbstractOptimizer alf::AlfonsoOpt + objsense::MOI.OptimizationSense + + function Optimizer(alf::AlfonsoOpt) + opt = new() + opt.alf = alf + return opt + end end Optimizer() = Optimizer(AlfonsoOpt()) @@ -263,13 +270,17 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ # TODO iterate over types of S, from easiest to hardest, appending cones in order of toughness # then SOC, RSOC, exp, power, sumofsquares, psd + + A = dropzeros!(sparse(IA, JA, VA, m, n)) b = Vector(dropzeros!(sparsevec(Ib, Vb, m))) c = Vector(dropzeros!(sparsevec(Jc, Vc, n))) - - # load problem data through native interface load_data!(opt.alf, A, b, c, cone) + opt.objsense = MOI.get(src, MOI.ObjectiveSense()) + + println(idxmap) + return idxmap end @@ -291,15 +302,32 @@ function MOI.get(opt::Optimizer, ::MOI.TerminationStatus) end end -MOI.get(opt::Optimizer, ::MOI.ObjectiveValue) = opt.alf.pobj -MOI.get(opt::Optimizer, ::MOI.ObjectiveBound) = opt.alf.dobj +function MOI.get(opt::Optimizer, ::MOI.ObjectiveValue) + if opt.objsense == MOI.MinSense + return get_pobj(opt.alf) + elseif opt.objsense == MOI.MaxSense + return -get_pobj(opt.alf) + else + return NaN + end +end + +function MOI.get(opt::Optimizer, ::MOI.ObjectiveBound) + if opt.objsense == MOI.MinSense + return get_dobj(opt.alf) + elseif opt.objsense == MOI.MaxSense + return -get_dobj(opt.alf) + else + return NaN + end +end function MOI.get(opt::Optimizer, ::MOI.ResultCount) - # TODO if status in (:Optimal, :NearlyInfeasible, :IllPosed) status = get_status(opt.alf) - if status == :Optimal + if status in (:Optimal, :NearlyInfeasible, :IllPosed) return 1 end + return 0 end function MOI.get(opt::Optimizer, ::MOI.PrimalStatus) @@ -319,3 +347,121 @@ function MOI.get(opt::Optimizer, ::MOI.DualStatus) return MOI.UnknownResultStatus end end + + +# TODO from ECOS: +# Implements getter for result value and statuses +# function MOI.get(instance::Optimizer, ::MOI.TerminationStatus) +# flag = instance.sol.ret_val +# if flag == ECOS.ECOS_OPTIMAL +# MOI.Success +# elseif flag == ECOS.ECOS_PINF +# MOI.Success +# elseif flag == ECOS.ECOS_DINF # Dual infeasible = primal unbounded, probably +# MOI.Success +# elseif flag == ECOS.ECOS_MAXIT +# MOI.IterationLimit +# elseif flag == ECOS.ECOS_OPTIMAL + ECOS.ECOS_INACC_OFFSET +# m.solve_stat = MOI.AlmostSuccess +# else +# m.solve_stat = MOI.OtherError +# end +# end +# +# MOI.get(instance::Optimizer, ::MOI.ObjectiveValue) = instance.sol.objval +# MOI.get(instance::Optimizer, ::MOI.ObjectiveBound) = instance.sol.objbnd +# +# function MOI.get(instance::Optimizer, ::MOI.PrimalStatus) +# flag = instance.sol.ret_val +# if flag == ECOS.ECOS_OPTIMAL +# MOI.FeasiblePoint +# elseif flag == ECOS.ECOS_PINF +# MOI.InfeasiblePoint +# elseif flag == ECOS.ECOS_DINF # Dual infeasible = primal unbounded, probably +# MOI.InfeasibilityCertificate +# elseif flag == ECOS.ECOS_MAXIT +# MOI.UnknownResultStatus +# elseif flag == ECOS.ECOS_OPTIMAL + ECOS.ECOS_INACC_OFFSET +# m.solve_stat = MOI.NearlyFeasiblePoint +# else +# m.solve_stat = MOI.OtherResultStatus +# end +# end +# function MOI.get(instance::Optimizer, ::MOI.DualStatus) +# flag = instance.sol.ret_val +# if flag == ECOS.ECOS_OPTIMAL +# MOI.FeasiblePoint +# elseif flag == ECOS.ECOS_PINF +# MOI.InfeasibilityCertificate +# elseif flag == ECOS.ECOS_DINF # Dual infeasible = primal unbounded, probably +# MOI.InfeasiblePoint +# elseif flag == ECOS.ECOS_MAXIT +# MOI.UnknownResultStatus +# elseif flag == ECOS.ECOS_OPTIMAL + ECOS.ECOS_INACC_OFFSET +# m.solve_stat = MOI.NearlyFeasiblePoint +# else +# m.solve_stat = MOI.OtherResultStatus +# end +# end + + +MOI.get(opt::Optimizer, ::MOI.VariablePrimal, vi::MOI.VariableIndex) = get_x(opt.alf)[vi.value] +MOI.get(opt::Optimizer, a::MOI.VariablePrimal, vi::Vector{MOI.VariableIndex}) = MOI.get.(opt, Ref(a), vi) + +# function MOI.get(opt::Optimizer, ::MOI.ConstraintPrimal, ci::MOI.ConstraintIndex{<:MOI.AbstractFunction, MOI.AbstractSet}) +# end +# MOI.get(opt::Optimizer, a::MOI.ConstraintPrimal, ci::Vector{MOI.ConstraintIndex}) = MOI.get.(opt, Ref(a), ci) + +MOI.get(opt::Optimizer, ::MOI.ConstraintDual, ci::MOI.ConstraintIndex) = get_y(opt.alf)[ci.value] +MOI.get(opt::Optimizer, a::MOI.ConstraintDual, ci::Vector{MOI.ConstraintIndex}) = MOI.get.(opt, Ref(a), ci) + + + + +# function MOI.get(instance::Optimizer, ::MOI.ConstraintPrimal, ci::CI{<:MOI.AbstractFunction, MOI.Zeros}) +# rows = constrrows(instance, ci) +# zeros(length(rows)) +# end +# function MOI.get(instance::Optimizer, ::MOI.ConstraintPrimal, ci::CI{<:MOI.AbstractFunction, <:MOI.EqualTo}) +# offset = constroffset(instance, ci) +# setconstant(instance, offset, ci) +# end +# function MOI.get(instance::Optimizer, ::MOI.ConstraintPrimal, ci::CI{<:MOI.AbstractFunction, S}) where S <: MOI.AbstractSet +# offset = constroffset(instance, ci) +# rows = constrrows(instance, ci) +# _unshift(instance, offset, scalecoef(rows, reorderval(instance.sol.slack[offset .+ rows], S), false, S), ci) +# end +# +# _dual(instance, ci::CI{<:MOI.AbstractFunction, <:ZeroCones}) = instance.sol.dual_eq +# _dual(instance, ci::CI) = instance.sol.dual_ineq +# function MOI.get(instance::Optimizer, ::MOI.ConstraintDual, ci::CI{<:MOI.AbstractFunction, S}) where S <: MOI.AbstractSet +# offset = constroffset(instance, ci) +# rows = constrrows(instance, ci) +# scalecoef(rows, reorderval(_dual(instance, ci)[offset .+ rows], S), false, S) +# end +# +# MOI.get(instance::Optimizer, ::MOI.ResultCount) = 1 + + + +# +# get_status(alf::AlfonsoOpt) = alf.status +# get_solvetime(alf::AlfonsoOpt) = alf.solvetime +# get_niters(alf::AlfonsoOpt) = alf.niters +# get_y(alf::AlfonsoOpt) = copy(alf.y) +# get_x(alf::AlfonsoOpt) = copy(alf.x) +# get_tau(alf::AlfonsoOpt) = alf.tau +# get_s(alf::AlfonsoOpt) = copy(alf.s) +# get_kappa(alf::AlfonsoOpt) = alf.kappa +# get_pobj(alf::AlfonsoOpt) = alf.pobj +# get_dobj(alf::AlfonsoOpt) = alf.dobj +# get_dgap(alf::AlfonsoOpt) = alf.dgap +# get_cgap(alf::AlfonsoOpt) = alf.cgap +# get_rel_dgap(alf::AlfonsoOpt) = alf.rel_dgap +# get_rel_cgap(alf::AlfonsoOpt) = alf.rel_cgap +# get_pres(alf::AlfonsoOpt) = copy(alf.pres) +# get_dres(alf::AlfonsoOpt) = copy(alf.dres) +# get_pin(alf::AlfonsoOpt) = alf.pin +# get_din(alf::AlfonsoOpt) = alf.din +# get_rel_pin(alf::AlfonsoOpt) = alf.rel_pin +# get_rel_din(alf::AlfonsoOpt) = alf.rel_din diff --git a/src/nativeinterface.jl b/src/nativeinterface.jl index e823bbcb7..15966060b 100644 --- a/src/nativeinterface.jl +++ b/src/nativeinterface.jl @@ -84,7 +84,7 @@ mutable struct AlfonsoOpt end function AlfonsoOpt(; - verbose = true, + verbose = false, optimtol = 1e-6, maxiter = 1e3, predlinesearch = true, @@ -116,10 +116,10 @@ end get_status(alf::AlfonsoOpt) = alf.status get_solvetime(alf::AlfonsoOpt) = alf.solvetime get_niters(alf::AlfonsoOpt) = alf.niters -get_y(alf::AlfonsoOpt) = copy(alf.y) -get_x(alf::AlfonsoOpt) = copy(alf.x) +get_y(alf::AlfonsoOpt) = alf.y # TODO? copy(alf.y) +get_x(alf::AlfonsoOpt) = alf.x # copy(alf.x) get_tau(alf::AlfonsoOpt) = alf.tau -get_s(alf::AlfonsoOpt) = copy(alf.s) +get_s(alf::AlfonsoOpt) = alf.s # copy(alf.s) get_kappa(alf::AlfonsoOpt) = alf.kappa get_pobj(alf::AlfonsoOpt) = alf.pobj get_dobj(alf::AlfonsoOpt) = alf.dobj @@ -127,8 +127,8 @@ get_dgap(alf::AlfonsoOpt) = alf.dgap get_cgap(alf::AlfonsoOpt) = alf.cgap get_rel_dgap(alf::AlfonsoOpt) = alf.rel_dgap get_rel_cgap(alf::AlfonsoOpt) = alf.rel_cgap -get_pres(alf::AlfonsoOpt) = copy(alf.pres) -get_dres(alf::AlfonsoOpt) = copy(alf.dres) +get_pres(alf::AlfonsoOpt) = alf.pres # copy(alf.pres) +get_dres(alf::AlfonsoOpt) = alf.dres # copy(alf.dres) get_pin(alf::AlfonsoOpt) = alf.pin get_din(alf::AlfonsoOpt) = alf.din get_rel_pin(alf::AlfonsoOpt) = alf.rel_pin @@ -208,7 +208,7 @@ function getinitialiterate(alf::AlfonsoOpt) cone = alf.cone # initial primal iterate - tx = similar(c) + tx = zeros(n) getintdir!(tx, cone) loadpnt!(cone, tx) @assert incone(cone) @@ -222,6 +222,7 @@ function getinitialiterate(alf::AlfonsoOpt) # central primal-dual iterate tx .*= sqrt(rp*rd) + @assert incone(cone) # TODO a little expensive to call this twice ty = zeros(m) tau = 1.0 diff --git a/test/runtests.jl b/test/runtests.jl index 2a6f44d95..d6154e934 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -38,9 +38,9 @@ const optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Opt const config = MOIT.TestConfig( atol=1e-4, rtol=1e-4, - solve=true, - query=false, - modify_lhs=false, + solve=false, + query=true, + modify_lhs=true, duals=false, infeas_certificates=false, ) @@ -60,3 +60,359 @@ end # )))), # config, exclude) # end + + +# +# +# +# model = optimizer +# +# atol = config.atol +# rtol = config.rtol +# # simple 2 variable, 1 constraint problem +# # min -x +# # st x + y <= 1 (x + y - 1 ∈ Nonpositives) +# # x, y >= 0 (x, y ∈ Nonnegatives) +# +# @test MOI.supports(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) +# @test MOI.supports(model, MOI.ObjectiveSense()) +# @test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}) +# @test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.GreaterThan{Float64}) +# @test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}) +# @test MOI.supports_constraint(model, MOI.SingleVariable, MOI.EqualTo{Float64}) +# @test MOI.supports_constraint(model, MOI.SingleVariable, MOI.GreaterThan{Float64}) +# +# #@test MOI.get(model, MOI.SupportsAddConstraintAfterSolve()) +# #@test MOI.get(model, MOI.SupportsAddVariableAfterSolve()) +# #@test MOI.get(model, MOI.SupportsDeleteConstraint()) +# +# MOI.empty!(model) +# @test MOI.is_empty(model) +# +# v = MOI.add_variables(model, 2) +# @test MOI.get(model, MOI.NumberOfVariables()) == 2 +# +# cf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0,1.0], v), 0.0) +# c = MOI.add_constraint(model, cf, MOI.LessThan(1.0)) +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}}()) == 1 +# +# vc1 = MOI.add_constraint(model, MOI.SingleVariable(v[1]), MOI.GreaterThan(0.0)) +# # test fallback +# vc2 = MOI.add_constraint(model, v[2], MOI.GreaterThan(0.0)) +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 2 +# +# # note: adding some redundant zero coefficients to catch solvers that don't handle duplicate coefficients correctly: +# objf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([0.0,0.0,-1.0,0.0,0.0,0.0], [v; v; v]), 0.0) +# MOI.set(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), objf) +# MOI.set(model, MOI.ObjectiveSense(), MOI.MinSense) +# +# @test MOI.get(model, MOI.ObjectiveSense()) == MOI.MinSense +# +# if config.query +# vrs = MOI.get(model, MOI.ListOfVariableIndices()) +# @test vrs == v || vrs == reverse(v) +# +# @test objf ≈ MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) +# +# @test cf ≈ MOI.get(model, MOI.ConstraintFunction(), c) +# +# s = MOI.get(model, MOI.ConstraintSet(), c) +# @test s == MOI.LessThan(1.0) +# +# s = MOI.get(model, MOI.ConstraintSet(), vc1) +# @test s == MOI.GreaterThan(0.0) +# +# s = MOI.get(model, MOI.ConstraintSet(), vc2) +# @test s == MOI.GreaterThan(0.0) +# end +# +# if config.solve +# MOI.optimize!(model) +# +# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success +# +# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint +# +# @test MOI.get(model, MOI.ObjectiveValue()) ≈ -1 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [1, 0] atol=atol rtol=rtol +# +# # @test MOI.get(model, MOI.ConstraintPrimal(), c) ≈ 1 atol=atol rtol=rtol +# +# if config.duals +# @test MOI.get(model, MOI.DualStatus()) == MOI.FeasiblePoint +# @test MOI.get(model, MOI.ConstraintDual(), c) ≈ -1 atol=atol rtol=rtol +# +# # reduced costs +# @test MOI.get(model, MOI.ConstraintDual(), vc1) ≈ 0 atol=atol rtol=rtol +# @test MOI.get(model, MOI.ConstraintDual(), vc2) ≈ 1 atol=atol rtol=rtol +# end +# end +# +# # change objective to Max +x +# +# objf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0,0.0], v), 0.0) +# MOI.set(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), objf) +# MOI.set(model, MOI.ObjectiveSense(), MOI.MaxSense) +# +# if config.query +# @test objf ≈ MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) +# end +# +# @test MOI.get(model, MOI.ObjectiveSense()) == MOI.MaxSense +# +# if config.solve +# MOI.optimize!(model) +# +# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success +# +# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint +# +# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 1 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [1, 0] atol=atol rtol=rtol +# +# if config.duals +# @test MOI.get(model, MOI.DualStatus()) == MOI.FeasiblePoint +# @test MOI.get(model, MOI.ConstraintDual(), c) ≈ -1 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.ConstraintDual(), vc1) ≈ 0 atol=atol rtol=rtol +# @test MOI.get(model, MOI.ConstraintDual(), vc2) ≈ 1 atol=atol rtol=rtol +# end +# end +# +# # add new variable to get : +# # max x + 2z +# # s.t. x + y + z <= 1 +# # x,y,z >= 0 +# +# z = MOI.add_variable(model) +# push!(v, z) +# @test v[3] == z +# +# if config.query +# # Test that the modification of v has not affected the model +# vars = map(t -> t.variable_index, MOI.get(model, MOI.ConstraintFunction(), c).terms) +# @test vars == [v[1], v[2]] || vars == [v[2], v[1]] +# @test MOI.ScalarAffineFunction([MOI.ScalarAffineTerm(1.0, v[1])], 0.0) ≈ MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) +# end +# +# vc3 = MOI.add_constraint(model, MOI.SingleVariable(v[3]), MOI.GreaterThan(0.0)) +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 3 +# +# if config.modify_lhs +# MOI.modify(model, c, MOI.ScalarCoefficientChange{Float64}(z, 1.0)) +# else +# MOI.delete(model, c) +# cf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0,1.0,1.0], v), 0.0) +# c = MOI.add_constraint(model, cf, MOI.LessThan(1.0)) +# end +# +# MOI.modify(model, +# MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), +# MOI.ScalarCoefficientChange{Float64}(z, 2.0) +# ) +# +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}}()) == 1 +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 3 +# +# if config.solve +# MOI.optimize!(model) +# +# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success +# +# @test MOI.get(model, MOI.ResultCount()) >= 1 +# +# @test MOI.get(model, MOI.PrimalStatus(1)) == MOI.FeasiblePoint +# +# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 2 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [0, 0, 1] atol=atol rtol=rtol +# +# # @test MOI.get(model, MOI.ConstraintPrimal(), c) ≈ 1 atol=atol rtol=rtol +# +# if config.duals +# @test MOI.get(model, MOI.DualStatus()) == MOI.FeasiblePoint +# @test MOI.get(model, MOI.ConstraintDual(), c) ≈ -2 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.ConstraintDual(), vc1) ≈ 1 atol=atol rtol=rtol +# @test MOI.get(model, MOI.ConstraintDual(), vc2) ≈ 2 atol=atol rtol=rtol +# @test MOI.get(model, MOI.ConstraintDual(), vc3) ≈ 0 atol=atol rtol=rtol +# end +# end +# +# # setting lb of x to -1 to get : +# # max x + 2z +# # s.t. x + y + z <= 1 +# # x >= -1 +# # y,z >= 0 +# MOI.set(model, MOI.ConstraintSet(), vc1, MOI.GreaterThan(-1.0)) +# +# if config.solve +# MOI.optimize!(model) +# +# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success +# +# @test MOI.get(model, MOI.ResultCount()) >= 1 +# +# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint +# +# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 3 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [-1, 0, 2] atol=atol rtol=rtol +# end +# +# # put lb of x back to 0 and fix z to zero to get : +# # max x + 2z +# # s.t. x + y + z <= 1 +# # x, y >= 0, z = 0 (vc3) +# MOI.set(model, MOI.ConstraintSet(), vc1, MOI.GreaterThan(0.0)) +# +# MOI.delete(model, vc3) +# +# vc3 = MOI.add_constraint(model, MOI.SingleVariable(v[3]), MOI.EqualTo(0.0)) +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 2 +# +# if config.solve +# MOI.optimize!(model) +# +# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success +# +# @test MOI.get(model, MOI.ResultCount()) >= 1 +# +# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint +# +# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 1 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [1, 0, 0] atol=atol rtol=rtol +# end +# +# # modify affine linear constraint set to be == 2 to get : +# # max x + 2z +# # s.t. x + y + z == 2 (c) +# # x,y >= 0, z = 0 +# MOI.delete(model, c) +# # note: adding some redundant zero coefficients to catch solvers that don't handle duplicate coefficients correctly: +# cf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([0.0,0.0,0.0,1.0,1.0,1.0,0.0,0.0,0.0], [v; v; v]), 0.0) +# c = MOI.add_constraint(model, cf, MOI.EqualTo(2.0)) +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}}()) == 0 +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.EqualTo{Float64}}()) == 1 +# +# if config.solve +# MOI.optimize!(model) +# +# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success +# +# @test MOI.get(model, MOI.ResultCount()) >= 1 +# +# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint +# +# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 2 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [2, 0, 0] atol=atol rtol=rtol +# end +# +# # modify objective function to x + 2y to get : +# # max x + 2y +# # s.t. x + y + z == 2 (c) +# # x,y >= 0, z = 0 +# +# objf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0,2.0,0.0], v), 0.0) +# MOI.set(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), objf) +# MOI.set(model, MOI.ObjectiveSense(), MOI.MaxSense) +# +# if config.query +# @test objf ≈ MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) +# end +# +# if config.solve +# MOI.optimize!(model) +# +# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success +# +# @test MOI.get(model, MOI.ResultCount()) >= 1 +# +# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint +# +# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 4 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [0, 2, 0] atol=atol rtol=rtol +# end +# +# # add constraint x - y >= 0 (c2) to get : +# # max x+2y +# # s.t. x + y + z == 2 (c) +# # x - y >= 0 (c2) +# # x,y >= 0 (vc1,vc2), z = 0 (vc3) +# +# cf2 = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0, -1.0, 0.0], v), 0.0) +# c2 = MOI.add_constraint(model, cf2, MOI.GreaterThan(0.0)) +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.EqualTo{Float64}}()) == 1 +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.GreaterThan{Float64}}()) == 1 +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}}()) == 0 +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.EqualTo{Float64}}()) == 1 +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 2 +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.LessThan{Float64}}()) == 0 +# +# if config.solve +# MOI.optimize!(model) +# +# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success +# +# @test MOI.get(model, MOI.ResultCount()) >= 1 +# +# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint +# +# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 3 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [1, 1, 0] atol=atol rtol=rtol +# +# # @test MOI.get(model, MOI.ConstraintPrimal(), c) ≈ 2 atol=atol rtol=rtol +# # +# # @test MOI.get(model, MOI.ConstraintPrimal(), c2) ≈ 0 atol=atol rtol=rtol +# # +# # @test MOI.get(model, MOI.ConstraintPrimal(), vc1) ≈ 1 atol=atol rtol=rtol +# # +# # @test MOI.get(model, MOI.ConstraintPrimal(), vc2) ≈ 1 atol=atol rtol=rtol +# # +# # @test MOI.get(model, MOI.ConstraintPrimal(), vc3) ≈ 0 atol=atol rtol=rtol +# +# if config.duals +# @test MOI.get(model, MOI.DualStatus(1)) == MOI.FeasiblePoint +# +# @test MOI.get(model, MOI.ConstraintDual(), c) ≈ -1.5 atol=atol rtol=rtol +# @test MOI.get(model, MOI.ConstraintDual(), c2) ≈ 0.5 atol=atol rtol=rtol +# +# @test MOI.get(model, MOI.ConstraintDual(), vc1) ≈ 0 atol=atol rtol=rtol +# @test MOI.get(model, MOI.ConstraintDual(), vc2) ≈ 0 atol=atol rtol=rtol +# @test MOI.get(model, MOI.ConstraintDual(), vc3) ≈ 1.5 atol=atol rtol=rtol +# end +# end +# +# if config.query +# @test MOI.get(model, MOI.ConstraintFunction(), c2) ≈ cf2 +# end +# +# # delete variable x to get : +# # max 2y +# # s.t. y + z == 2 +# # - y >= 0 +# # y >= 0, z = 0 +# +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 2 +# MOI.delete(model, v[1]) +# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 1 +# +# if config.query +# @test MOI.get(model, MOI.ConstraintFunction(), c2) ≈ MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([-1.0, 0.0], [v[2], z]), 0.0) +# +# vrs = MOI.get(model, MOI.ListOfVariableIndices()) +# @test vrs == [v[2], z] || vrs == [z, v[2]] +# @test MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) ≈ MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([2.0, 0.0], [v[2], z]), 0.0) +# end +# +# + + + +return nothing From f0521122f60d5e60e328a8df236f4f64f0939829 Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Tue, 11 Sep 2018 17:37:16 -0400 Subject: [PATCH 04/10] support MOI cones --- src/mathoptinterface.jl | 192 ++++++++++++++++++++++++---------------- test/runtests.jl | 30 +++---- 2 files changed, 129 insertions(+), 93 deletions(-) diff --git a/src/mathoptinterface.jl b/src/mathoptinterface.jl index 3c45f3e83..073b99779 100644 --- a/src/mathoptinterface.jl +++ b/src/mathoptinterface.jl @@ -44,14 +44,6 @@ const SupportedSets = Union{ MOI.supports_constraint(::Optimizer, ::Type{<:SupportedFuns}, ::Type{<:SupportedSets}) = true -conefrommoi(s::MOI.Nonnegatives) = NonnegativeCone(MOI.dimension(s)) -conefrommoi(s::MOI.Nonpositives) = error("Nonpositive cone not handled yet") -conefrommoi(s::MOI.SecondOrderCone) = SecondOrderCone(MOI.dimension(s)) -conefrommoi(s::MOI.RotatedSecondOrderCone) = RotatedSecondOrderCone(MOI.dimension(s)) -conefrommoi(s::MOI.PositiveSemidefiniteConeTriangle) = PositiveSemidefiniteCone(MOI.dimension(s)) -conefrommoi(s::MOI.ExponentialCone) = ExponentialCone() -conefrommoi(s::MOI.PowerCone) = PowerCone(s.exponent) - # build representation as min c'x s.t. Ax = b, x in K # TODO what if some variables x are in multiple cone constraints? function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_attributes=true) @@ -183,23 +175,46 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ nonnegvars = Int[] # for building up a single nonnegative cone - for S in (MOI.GreaterThan{Float64}, MOI.LessThan{Float64}), ci in getsrccons(MOI.SingleVariable, S) - i += 1 - idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, S}(i) - - bval = ((S == MOI.GreaterThan{Float64}) ? getconset(ci).lower : getconset(ci).upper) - vn = getconfun(ci).variable - if iszero(bval) - # no need to change A,b - push!(nonnegvars, idxmap[vn].value) - else + for S in (MOI.GreaterThan{Float64}, MOI.LessThan{Float64}) + for ci in getsrccons(MOI.SingleVariable, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, S}(i) + + bval = ((S == MOI.GreaterThan{Float64}) ? getconset(ci).lower : getconset(ci).upper) + vn = getconfun(ci).variable + if iszero(bval) + # no need to change A,b + push!(nonnegvars, idxmap[vn].value) + else + # add auxiliary variable and equality constraint + m += 1 + push!(Ib, m) + push!(Vb, bval) + push!(IA, m) + push!(JA, idxmap[vn].value) + push!(VA, 1.0) + n += 1 + push!(IA, m) + push!(JA, n) + push!(VA, ((S == MOI.GreaterThan{Float64}) ? -1.0 : 1.0)) + push!(nonnegvars, n) + end + end + + for ci in getsrccons(MOI.ScalarAffineFunction{Float64}, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, S}(i) + # add auxiliary variable and equality constraint + fi = getconfun(ci) m += 1 push!(Ib, m) - push!(Vb, bval) - push!(IA, m) - push!(JA, idxmap[vn].value) - push!(VA, 1.0) + push!(Vb, ((S == MOI.GreaterThan{Float64}) ? getconset(ci).lower : getconset(ci).upper) - fi.constant) + for vt in fi.terms + push!(IA, m) + push!(JA, idxmap[vt.variable_index].value) + push!(VA, vt.coefficient) + end n += 1 push!(IA, m) push!(JA, n) @@ -208,79 +223,100 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ end end - for S in (MOI.GreaterThan{Float64}, MOI.LessThan{Float64}), ci in getsrccons(MOI.ScalarAffineFunction{Float64}, S) - i += 1 - idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, S}(i) - - # add auxiliary variable and equality constraint - fi = getconfun(ci) - m += 1 - push!(Ib, m) - push!(Vb, ((S == MOI.GreaterThan{Float64}) ? getconset(ci).lower : getconset(ci).upper) - fi.constant) - for vt in fi.terms - push!(IA, m) - push!(JA, idxmap[vt.variable_index].value) - push!(VA, vt.coefficient) - end - n += 1 - push!(IA, m) - push!(JA, n) - push!(VA, ((S == MOI.GreaterThan{Float64}) ? -1.0 : 1.0)) - push!(nonnegvars, n) - end - - for S in (MOI.Nonnegatives, MOI.Nonpositives), ci in getsrccons(MOI.VectorOfVariables, S) - i += 1 - idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, S}(i) - - for vj in getconfun(ci).variables - col = idxmap[vj].value - m += 1 - push!(IA, m) - push!(JA, col) - push!(VA, ((S == MOI.Nonnegatives) ? 1.0 : -1.0)) - push!(nonnegvars, col) + for S in (MOI.Nonnegatives, MOI.Nonpositives) + for ci in getsrccons(MOI.VectorOfVariables, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, S}(i) + + for vj in getconfun(ci).variables + col = idxmap[vj].value + m += 1 + push!(IA, m) + push!(JA, col) + push!(VA, ((S == MOI.Nonnegatives) ? 1.0 : -1.0)) + push!(nonnegvars, col) + end end - end - - for S in (MOI.Nonnegatives, MOI.Nonpositives), ci in getsrccons(MOI.VectorAffineFunction{Float64}, S) - i += 1 - idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, S}(i) - fi = getconfun(ci) - dim = MOI.output_dimension(fi) - for vt in fi.terms - push!(IA, m + vt.output_index) - push!(JA, idxmap[vt.scalar_term.variable_index].value) - push!(VA, vt.scalar_term.coefficient) + for ci in getsrccons(MOI.VectorAffineFunction{Float64}, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, S}(i) + + fi = getconfun(ci) + dim = MOI.output_dimension(fi) + for vt in fi.terms + push!(IA, m + vt.output_index) + push!(JA, idxmap[vt.scalar_term.variable_index].value) + push!(VA, vt.scalar_term.coefficient) + end + append!(Ib, collect(m+1:m+dim)) + append!(Vb, -fi.constants) + append!(IA, collect(m+1:m+dim)) + append!(JA, collect(n+1:n+dim)) + append!(VA, fill(((S == MOI.Nonnegatives) ? 1.0 : -1.0), dim)) + m += dim + append!(nonnegvars, collect(n+1:n+dim)) + n += dim end - append!(Ib, collect(m+1:m+dim)) - append!(Vb, -fi.constants) - append!(IA, collect(m+1:m+dim)) - append!(JA, collect(n+1:n+dim)) - append!(VA, fill(((S == MOI.Nonnegatives) ? 1.0 : -1.0), dim)) - m += dim - n += dim end # add single nonnegative cone addprimitivecone!(cone, NonnegativeCone(length(nonnegvars)), nonnegvars) + # add other primitive cone constraints (in increasing order of difficulty of evaluating in-cone checks) + + conefrommoi(s::MOI.SecondOrderCone) = SecondOrderCone(MOI.dimension(s)) + conefrommoi(s::MOI.RotatedSecondOrderCone) = RotatedSecondOrderCone(MOI.dimension(s)) + conefrommoi(s::MOI.PositiveSemidefiniteConeTriangle) = PositiveSemidefiniteCone(MOI.dimension(s)) + conefrommoi(s::MOI.ExponentialCone) = ExponentialCone() + conefrommoi(s::MOI.PowerCone) = PowerCone(s.exponent) + + for S in (MOI.SecondOrderCone, MOI.RotatedSecondOrderCone, MOI.ExponentialCone, MOI.PowerCone, MOI.PositiveSemidefiniteConeTriangle) + for ci in getsrccons(MOI.VectorOfVariables, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, S}(i) + + fi = getconfun(ci) + for vj in fi.variables + m += 1 + push!(IA, m) + push!(JA, idxmap[vj].value) + push!(VA, 1.0) + end + addprimitivecone!(cone, conefrommoi(getconset(ci)), [idxmap[vj].value for vj in fi.variables]) + end - # TODO iterate over types of S, from easiest to hardest, appending cones in order of toughness - # then SOC, RSOC, exp, power, sumofsquares, psd - - + for ci in getsrccons(MOI.VectorAffineFunction{Float64}, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, S}(i) + + fi = getconfun(ci) + dim = MOI.output_dimension(fi) + for vt in fi.terms + push!(IA, m + vt.output_index) + push!(JA, idxmap[vt.scalar_term.variable_index].value) + push!(VA, vt.scalar_term.coefficient) + end + append!(Ib, collect(m+1:m+dim)) + append!(Vb, -fi.constants) + append!(IA, collect(m+1:m+dim)) + append!(JA, collect(n+1:n+dim)) + append!(VA, fill(1.0, dim)) + m += dim + addprimitivecone!(cone, conefrommoi(getconset(ci)), n+1:n+dim) + n += dim + end + end + # finalize data and load into Alfonso model A = dropzeros!(sparse(IA, JA, VA, m, n)) b = Vector(dropzeros!(sparsevec(Ib, Vb, m))) c = Vector(dropzeros!(sparsevec(Jc, Vc, n))) load_data!(opt.alf, A, b, c, cone) + # set objective sense opt.objsense = MOI.get(src, MOI.ObjectiveSense()) - println(idxmap) - return idxmap end diff --git a/test/runtests.jl b/test/runtests.jl index d6154e934..d063c14d9 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -38,29 +38,29 @@ const optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Opt const config = MOIT.TestConfig( atol=1e-4, rtol=1e-4, - solve=false, + solve=true, query=true, modify_lhs=true, duals=false, infeas_certificates=false, ) -@testset "Continuous linear problems" begin - MOIT.contlineartest(MOIB.SplitInterval{Float64}(optimizer), config) -end - -# @testset "Continuous conic problems" begin -# exclude = ["rootdet", "logdet"] -# MOIT.contconictest( -# MOIB.SquarePSD{Float64}( -# MOIB.GeoMean{Float64}( -# MOIB.LogDet{Float64}( -# MOIB.RootDet{Float64}( -# optimizer -# )))), -# config, exclude) +# @testset "Continuous linear problems" begin +# MOIT.contlineartest(MOIB.SplitInterval{Float64}(optimizer), config) # end +@testset "Continuous conic problems" begin + exclude = ["rootdet", "logdet"] + MOIT.contconictest( + MOIB.SquarePSD{Float64}( + MOIB.GeoMean{Float64}( + # MOIB.LogDet{Float64}( + # MOIB.RootDet{Float64}( + optimizer + )), + config, exclude) +end + # # From 35e042ff41a9d576fdd5c53ec76c01e45367ee5f Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Sat, 29 Sep 2018 14:56:05 -0400 Subject: [PATCH 05/10] add nonpositives cone, start updating MOI functions for new P-D format --- src/Alfonso.jl | 8 +-- src/cone.jl | 19 ------ src/mathoptinterface.jl | 94 ++++++++++++++------------ src/primitivecones/nonnegative.jl | 31 --------- src/primitivecones/orthant.jl | 53 +++++++++++++++ test/nativeexamples.jl | 107 ++++++++++++++++-------------- test/runtests.jl | 106 ++++++++++++++--------------- 7 files changed, 219 insertions(+), 199 deletions(-) delete mode 100644 src/primitivecones/nonnegative.jl create mode 100644 src/primitivecones/orthant.jl diff --git a/src/Alfonso.jl b/src/Alfonso.jl index c8948cebe..403b4f936 100644 --- a/src/Alfonso.jl +++ b/src/Alfonso.jl @@ -7,7 +7,7 @@ module Alfonso include("interpolation.jl") include("cone.jl") for primcone in [ - "nonnegative", + "orthant", "sumofsquares", "secondorder", "exponential", @@ -20,7 +20,7 @@ module Alfonso end include("nativeinterface.jl") - import MathOptInterface - const MOI = MathOptInterface - include("mathoptinterface.jl") + # import MathOptInterface + # const MOI = MathOptInterface + # include("mathoptinterface.jl") end diff --git a/src/cone.jl b/src/cone.jl index 295c09631..8d941d266 100644 --- a/src/cone.jl +++ b/src/cone.jl @@ -9,25 +9,6 @@ abstract type PrimitiveCone end mutable struct Cone prms::Vector{PrimitiveCone} idxs::Vector{AbstractVector{Int}} - - function Cone(prms::Vector{PrimitiveCone}, idxs::Vector{AbstractVector{Int}}) - @assert length(prms) == length(idxs) - for k in eachindex(prms) - @assert dimension(prms[k]) == length(idxs[k]) - end - cone = new() - cone.prms = prms - cone.idxs = idxs - return cone - end -end -Cone() = Cone(PrimitiveCone[], AbstractVector{Int}[]) - -function addprimitivecone!(cone::Cone, prm::PrimitiveCone, idx::AbstractVector{Int}) - @assert dimension(prm) == length(idx) - push!(cone.prms, prm) - push!(cone.idxs, idx) - return cone end # calculate complexity parameter of the barrier (sum of the primitive cone barrier parameters) diff --git a/src/mathoptinterface.jl b/src/mathoptinterface.jl index 073b99779..8ba97f1d3 100644 --- a/src/mathoptinterface.jl +++ b/src/mathoptinterface.jl @@ -39,7 +39,7 @@ const SupportedSets = Union{ MOI.RotatedSecondOrderCone, MOI.PositiveSemidefiniteConeTriangle, MOI.ExponentialCone, - MOI.PowerCone, + # MOI.PowerCone, } MOI.supports_constraint(::Optimizer, ::Type{<:SupportedFuns}, ::Type{<:SupportedSets}) = true @@ -79,36 +79,34 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ getconset(conidx) = MOI.get(src, MOI.ConstraintSet(), conidx) # pass over vector constraints to construct conic model data - # TODO don't enforce Float64 type (IA, JA, VA) = (Int[], Int[], Float64[]) (Ib, Vb) = (Int[], Float64[]) cone = Cone() i = 0 # MOI constraint objects - m = 0 # rows of A + p = 0 # rows of A (equality constraint matrix) + q = 0 # rows of G (cone constraint matrix) # equality constraints for ci in getsrccons(MOI.SingleVariable, MOI.EqualTo{Float64}) - # TODO can preprocess out maybe i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, MOI.EqualTo{Float64}}(i) - m += 1 - push!(Ib, m) + p += 1 + push!(Ib, p) push!(Vb, getconset(ci).value) - push!(IA, m) + push!(IA, p) push!(JA, idxmap[getconfun(ci).variable].value) push!(VA, 1.0) end for ci in getsrccons(MOI.SingleVariable, MOI.Zeros) - # TODO can preprocess out maybe i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, MOI.Zeros}(i) - m += 1 - push!(IA, m) + p += 1 + push!(IA, p) push!(JA, idxmap[getconfun(ci).variable].value) push!(VA, 1.0) end @@ -118,11 +116,11 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}}(i) fi = getconfun(ci) - m += 1 - push!(Ib, m) + p += 1 + push!(Ib, p) push!(Vb, getconset(ci).value - fi.constant) for vt in fi.terms - push!(IA, m) + push!(IA, p) push!(JA, idxmap[vt.variable_index].value) push!(VA, vt.coefficient) end @@ -133,11 +131,11 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, MOI.Zeros}(i) fi = getconfun(ci) - m += 1 - push!(Ib, m) + p += 1 + push!(Ib, p) push!(Vb, -fi.constant) for vt in fi.terms - push!(IA, m) + push!(IA, p) push!(JA, idxmap[vt.variable_index].value) push!(VA, vt.coefficient) end @@ -148,8 +146,8 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, MOI.Zeros}(i) for vi in getconfun(ci).variables - m += 1 - push!(IA, m) + p += 1 + push!(IA, p) push!(JA, idxmap[vi].value) push!(VA, 1.0) end @@ -161,19 +159,27 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ fi = getconfun(ci) dim = MOI.output_dimension(fi) - append!(Ib, collect(m+1:m+dim)) + append!(Ib, collect(p+1:p+dim)) append!(Vb, -fi.constants) for vt in fi.terms - push!(IA, m + vt.output_index) + push!(IA, p + vt.output_index) push!(JA, idxmap[vt.scalar_term.variable_index].value) push!(VA, vt.scalar_term.coefficient) end - m += dim + p += dim end + + + + + + + # linear inequality constraints nonnegvars = Int[] # for building up a single nonnegative cone + # TODO do the same for variable bounds to make a L_inf cone for S in (MOI.GreaterThan{Float64}, MOI.LessThan{Float64}) for ci in getsrccons(MOI.SingleVariable, S) @@ -187,14 +193,14 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ push!(nonnegvars, idxmap[vn].value) else # add auxiliary variable and equality constraint - m += 1 - push!(Ib, m) + p += 1 + push!(Ib, p) push!(Vb, bval) - push!(IA, m) + push!(IA, p) push!(JA, idxmap[vn].value) push!(VA, 1.0) n += 1 - push!(IA, m) + push!(IA, p) push!(JA, n) push!(VA, ((S == MOI.GreaterThan{Float64}) ? -1.0 : 1.0)) push!(nonnegvars, n) @@ -207,16 +213,16 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ # add auxiliary variable and equality constraint fi = getconfun(ci) - m += 1 - push!(Ib, m) + p += 1 + push!(Ib, p) push!(Vb, ((S == MOI.GreaterThan{Float64}) ? getconset(ci).lower : getconset(ci).upper) - fi.constant) for vt in fi.terms - push!(IA, m) + push!(IA, p) push!(JA, idxmap[vt.variable_index].value) push!(VA, vt.coefficient) end n += 1 - push!(IA, m) + push!(IA, p) push!(JA, n) push!(VA, ((S == MOI.GreaterThan{Float64}) ? -1.0 : 1.0)) push!(nonnegvars, n) @@ -230,8 +236,8 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for vj in getconfun(ci).variables col = idxmap[vj].value - m += 1 - push!(IA, m) + p += 1 + push!(IA, p) push!(JA, col) push!(VA, ((S == MOI.Nonnegatives) ? 1.0 : -1.0)) push!(nonnegvars, col) @@ -245,16 +251,16 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ fi = getconfun(ci) dim = MOI.output_dimension(fi) for vt in fi.terms - push!(IA, m + vt.output_index) + push!(IA, p + vt.output_index) push!(JA, idxmap[vt.scalar_term.variable_index].value) push!(VA, vt.scalar_term.coefficient) end - append!(Ib, collect(m+1:m+dim)) + append!(Ib, collect(p+1:p+dim)) append!(Vb, -fi.constants) - append!(IA, collect(m+1:m+dim)) + append!(IA, collect(p+1:p+dim)) append!(JA, collect(n+1:n+dim)) append!(VA, fill(((S == MOI.Nonnegatives) ? 1.0 : -1.0), dim)) - m += dim + p += dim append!(nonnegvars, collect(n+1:n+dim)) n += dim end @@ -269,7 +275,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ conefrommoi(s::MOI.RotatedSecondOrderCone) = RotatedSecondOrderCone(MOI.dimension(s)) conefrommoi(s::MOI.PositiveSemidefiniteConeTriangle) = PositiveSemidefiniteCone(MOI.dimension(s)) conefrommoi(s::MOI.ExponentialCone) = ExponentialCone() - conefrommoi(s::MOI.PowerCone) = PowerCone(s.exponent) + # conefrommoi(s::MOI.PowerCone) = PowerCone(s.exponent) for S in (MOI.SecondOrderCone, MOI.RotatedSecondOrderCone, MOI.ExponentialCone, MOI.PowerCone, MOI.PositiveSemidefiniteConeTriangle) for ci in getsrccons(MOI.VectorOfVariables, S) @@ -278,8 +284,8 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ fi = getconfun(ci) for vj in fi.variables - m += 1 - push!(IA, m) + p += 1 + push!(IA, p) push!(JA, idxmap[vj].value) push!(VA, 1.0) end @@ -293,24 +299,24 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ fi = getconfun(ci) dim = MOI.output_dimension(fi) for vt in fi.terms - push!(IA, m + vt.output_index) + push!(IA, p + vt.output_index) push!(JA, idxmap[vt.scalar_term.variable_index].value) push!(VA, vt.scalar_term.coefficient) end - append!(Ib, collect(m+1:m+dim)) + append!(Ib, collect(p+1:p+dim)) append!(Vb, -fi.constants) - append!(IA, collect(m+1:m+dim)) + append!(IA, collect(p+1:p+dim)) append!(JA, collect(n+1:n+dim)) append!(VA, fill(1.0, dim)) - m += dim + p += dim addprimitivecone!(cone, conefrommoi(getconset(ci)), n+1:n+dim) n += dim end end # finalize data and load into Alfonso model - A = dropzeros!(sparse(IA, JA, VA, m, n)) - b = Vector(dropzeros!(sparsevec(Ib, Vb, m))) + A = dropzeros!(sparse(IA, JA, VA, p, n)) + b = Vector(dropzeros!(sparsevec(Ib, Vb, p))) c = Vector(dropzeros!(sparsevec(Jc, Vc, n))) load_data!(opt.alf, A, b, c, cone) diff --git a/src/primitivecones/nonnegative.jl b/src/primitivecones/nonnegative.jl deleted file mode 100644 index 6f0e5edd0..000000000 --- a/src/primitivecones/nonnegative.jl +++ /dev/null @@ -1,31 +0,0 @@ - -# nonnegative orthant cone -# barrier is -sum_j ln x_j -# from Nesterov & Todd "Self-Scaled Barriers and Interior-Point Methods for Convex Programming" -mutable struct NonnegativeCone <: PrimitiveCone - dim::Int - pnt::AbstractVector{Float64} - invpnt::Vector{Float64} - - function NonnegativeCone(dim::Int) - prm = new() - prm.dim = dim - prm.invpnt = Vector{Float64}(undef, dim) - return prm - end -end - -dimension(prm::NonnegativeCone) = prm.dim -barrierpar_prm(prm::NonnegativeCone) = prm.dim -getintdir_prm!(arr::AbstractVector{Float64}, prm::NonnegativeCone) = (@. arr = 1.0; arr) -loadpnt_prm!(prm::NonnegativeCone, pnt::AbstractVector{Float64}) = (prm.pnt = pnt) -incone_prm(prm::NonnegativeCone) = all(x -> (x > 0.0), prm.pnt) - -function calcg_prm!(g::AbstractVector{Float64}, prm::NonnegativeCone) - @. prm.invpnt = inv(prm.pnt) - @. g = -prm.invpnt - return g -end - -calcHiarr_prm!(prod::AbstractArray{Float64}, arr::AbstractArray{Float64}, prm::NonnegativeCone) = (@. prod = abs2(prm.pnt)*arr; prod) -calcHarr_prm!(prod::AbstractArray{Float64}, arr::AbstractArray{Float64}, prm::NonnegativeCone) = (@. prod = abs2(prm.invpnt)*arr; prod) diff --git a/src/primitivecones/orthant.jl b/src/primitivecones/orthant.jl new file mode 100644 index 000000000..ba370fa83 --- /dev/null +++ b/src/primitivecones/orthant.jl @@ -0,0 +1,53 @@ + +# nonnegative/nonpositive orthant cones +# from Nesterov & Todd "Self-Scaled Barriers and Interior-Point Methods for Convex Programming" + +# nonnegative cone barrier is -sum_j ln x_j +mutable struct NonnegativeCone <: PrimitiveCone + dim::Int + pnt::AbstractVector{Float64} + invpnt::Vector{Float64} + + function NonnegativeCone(dim::Int) + prm = new() + prm.dim = dim + prm.invpnt = Vector{Float64}(undef, dim) + return prm + end +end + +# nonpositive cone barrier is -sum_j ln x_j +mutable struct NonpositiveCone <: PrimitiveCone + dim::Int + pnt::AbstractVector{Float64} + invpnt::Vector{Float64} + + function NonpositiveCone(dim::Int) + prm = new() + prm.dim = dim + prm.invpnt = Vector{Float64}(undef, dim) + return prm + end +end + +OrthantCone = Union{NonnegativeCone, NonpositiveCone} + +dimension(prm::OrthantCone) = prm.dim +barrierpar_prm(prm::OrthantCone) = prm.dim + +getintdir_prm!(arr::AbstractVector{Float64}, prm::NonnegativeCone) = (@. arr = 1.0; arr) +getintdir_prm!(arr::AbstractVector{Float64}, prm::NonpositiveCone) = (@. arr = -1.0; arr) + +loadpnt_prm!(prm::OrthantCone, pnt::AbstractVector{Float64}) = (prm.pnt = pnt) + +incone_prm(prm::NonnegativeCone) = all(x -> (x > 0.0), prm.pnt) +incone_prm(prm::NonpositiveCone) = all(x -> (x < 0.0), prm.pnt) + +function calcg_prm!(g::AbstractVector{Float64}, prm::OrthantCone) + @. prm.invpnt = inv(prm.pnt) + @. g = -prm.invpnt + return g +end + +calcHiarr_prm!(prod::AbstractArray{Float64}, arr::AbstractArray{Float64}, prm::OrthantCone) = (@. prod = abs2(prm.pnt)*arr; prod) +calcHarr_prm!(prod::AbstractArray{Float64}, arr::AbstractArray{Float64}, prm::OrthantCone) = (@. prod = abs2(prm.invpnt)*arr; prod) diff --git a/test/nativeexamples.jl b/test/nativeexamples.jl index afad7af93..877fec599 100644 --- a/test/nativeexamples.jl +++ b/test/nativeexamples.jl @@ -1,17 +1,28 @@ -@testset "small lp 1" begin - alf = Alfonso.AlfonsoOpt(verbose=verbflag) +@testset "small lp 1: nonnegative vs nonpositive orthant" begin Random.seed!(1) (n, p, q) = (10, 8, 10) c = rand(0.0:9.0, n) A = rand(-9.0:9.0, p, n) b = A*ones(n) - G = SparseMatrixCSC(-1.0I, q, n) h = zeros(q) + + alf1 = Alfonso.AlfonsoOpt(verbose=verbflag) + G = SparseMatrixCSC(-1.0I, q, n) cone = Alfonso.Cone([Alfonso.NonnegativeCone(q)], [1:q]) - Alfonso.load_data!(alf, c, A, b, G, h, cone) - @time Alfonso.solve!(alf) - @test Alfonso.get_status(alf) == :Optimal + Alfonso.load_data!(alf1, c, A, b, G, h, cone) + @time Alfonso.solve!(alf1) + @test Alfonso.get_status(alf1) == :Optimal + + alf2 = Alfonso.AlfonsoOpt(verbose=verbflag) + G = SparseMatrixCSC(1.0I, q, n) + cone = Alfonso.Cone([Alfonso.NonpositiveCone(q)], [1:q]) + Alfonso.load_data!(alf2, c, A, b, G, h, cone) + @time Alfonso.solve!(alf2) + @test Alfonso.get_status(alf2) == :Optimal + + @test Alfonso.get_pobj(alf1) ≈ Alfonso.get_pobj(alf2) atol=1e-4 rtol=1e-4 + @test Alfonso.get_dobj(alf1) ≈ Alfonso.get_dobj(alf2) atol=1e-4 rtol=1e-4 end @testset "small lp 2" begin @@ -37,9 +48,9 @@ end A = rand(-9.0:9.0, p, n) b = A*ones(n) @assert n == q - G = Diagonal(-1.0I, n) + G = Diagonal(1.0I, n) h = zeros(q) - cone = Alfonso.Cone([Alfonso.NonnegativeCone(q)], [1:q]) + cone = Alfonso.Cone([Alfonso.NonpositiveCone(q)], [1:q]) Alfonso.load_data!(alf, c, A, b, G, h, cone) @time Alfonso.solve!(alf) @test Alfonso.get_status(alf) == :Optimal @@ -202,41 +213,41 @@ end @testset "small dense lp example (dense vs sparse A)" begin # dense methods - d_alf = Alfonso.AlfonsoOpt(verbose=verbflag) - build_lp!(d_alf, 50, 100, dense=true, tosparse=false) - @time Alfonso.solve!(d_alf) - @test Alfonso.get_niters(d_alf) <= 40 - @test Alfonso.get_status(d_alf) == :Optimal + alf2 = Alfonso.AlfonsoOpt(verbose=verbflag) + build_lp!(alf2, 50, 100, dense=true, tosparse=false) + @time Alfonso.solve!(alf2) + @test Alfonso.get_niters(alf2) <= 40 + @test Alfonso.get_status(alf2) == :Optimal # sparse methods - s_alf = Alfonso.AlfonsoOpt(verbose=verbflag) - build_lp!(s_alf, 50, 100, dense=true, tosparse=true) - @time Alfonso.solve!(s_alf) - @test Alfonso.get_niters(s_alf) <= 40 - @test Alfonso.get_status(s_alf) == :Optimal - - @test Alfonso.get_pobj(d_alf) ≈ Alfonso.get_pobj(s_alf) atol=1e-4 rtol=1e-4 - @test Alfonso.get_dobj(d_alf) ≈ Alfonso.get_dobj(s_alf) atol=1e-4 rtol=1e-4 + alf1 = Alfonso.AlfonsoOpt(verbose=verbflag) + build_lp!(alf1, 50, 100, dense=true, tosparse=true) + @time Alfonso.solve!(alf1) + @test Alfonso.get_niters(alf1) <= 40 + @test Alfonso.get_status(alf1) == :Optimal + + @test Alfonso.get_pobj(alf2) ≈ Alfonso.get_pobj(alf1) atol=1e-4 rtol=1e-4 + @test Alfonso.get_dobj(alf2) ≈ Alfonso.get_dobj(alf1) atol=1e-4 rtol=1e-4 end @testset "1D poly envelope example (dense vs sparse A)" begin # dense methods - d_alf = Alfonso.AlfonsoOpt(verbose=verbflag) - build_envelope!(d_alf, 2, 5, 1, 5, use_data=true, dense=true) - @time Alfonso.solve!(d_alf) - @test Alfonso.get_niters(d_alf) <= 30 - @test Alfonso.get_status(d_alf) == :Optimal - @test Alfonso.get_pobj(d_alf) ≈ -25.502777 atol=1e-4 rtol=1e-4 - @test Alfonso.get_dobj(d_alf) ≈ -25.502777 atol=1e-4 rtol=1e-4 + alf2 = Alfonso.AlfonsoOpt(verbose=verbflag) + build_envelope!(alf2, 2, 5, 1, 5, use_data=true, dense=true) + @time Alfonso.solve!(alf2) + @test Alfonso.get_niters(alf2) <= 30 + @test Alfonso.get_status(alf2) == :Optimal + @test Alfonso.get_pobj(alf2) ≈ -25.502777 atol=1e-4 rtol=1e-4 + @test Alfonso.get_dobj(alf2) ≈ -25.502777 atol=1e-4 rtol=1e-4 # sparse methods - s_alf = Alfonso.AlfonsoOpt(verbose=verbflag) - build_envelope!(s_alf, 2, 5, 1, 5, use_data=true, dense=false) - @time Alfonso.solve!(s_alf) - @test Alfonso.get_niters(s_alf) <= 30 - @test Alfonso.get_status(s_alf) == :Optimal - @test Alfonso.get_pobj(s_alf) ≈ -25.502777 atol=1e-4 rtol=1e-4 - @test Alfonso.get_dobj(s_alf) ≈ -25.502777 atol=1e-4 rtol=1e-4 + alf1 = Alfonso.AlfonsoOpt(verbose=verbflag) + build_envelope!(alf1, 2, 5, 1, 5, use_data=true, dense=false) + @time Alfonso.solve!(alf1) + @test Alfonso.get_niters(alf1) <= 30 + @test Alfonso.get_status(alf1) == :Optimal + @test Alfonso.get_pobj(alf1) ≈ -25.502777 atol=1e-4 rtol=1e-4 + @test Alfonso.get_dobj(alf1) ≈ -25.502777 atol=1e-4 rtol=1e-4 end # most values taken from https://people.sc.fsu.edu/~jburkardt/py_src/polynomials/polynomials.html @@ -371,21 +382,21 @@ end @testset "2D poly envelope example (dense vs sparse A)" begin # dense methods - d_alf = Alfonso.AlfonsoOpt(verbose=verbflag) - build_envelope!(d_alf, 2, 4, 2, 7, dense=true) - @time Alfonso.solve!(d_alf) - @test Alfonso.get_niters(d_alf) <= 55 - @test Alfonso.get_status(d_alf) == :Optimal + alf2 = Alfonso.AlfonsoOpt(verbose=verbflag) + build_envelope!(alf2, 2, 4, 2, 7, dense=true) + @time Alfonso.solve!(alf2) + @test Alfonso.get_niters(alf2) <= 55 + @test Alfonso.get_status(alf2) == :Optimal # sparse methods - s_alf = Alfonso.AlfonsoOpt(verbose=verbflag) - build_envelope!(s_alf, 2, 4, 2, 7, dense=false) - @time Alfonso.solve!(s_alf) - @test Alfonso.get_niters(s_alf) <= 55 - @test Alfonso.get_status(s_alf) == :Optimal - - @test Alfonso.get_pobj(d_alf) ≈ Alfonso.get_pobj(s_alf) atol=1e-4 rtol=1e-4 - @test Alfonso.get_dobj(d_alf) ≈ Alfonso.get_dobj(s_alf) atol=1e-4 rtol=1e-4 + alf1 = Alfonso.AlfonsoOpt(verbose=verbflag) + build_envelope!(alf1, 2, 4, 2, 7, dense=false) + @time Alfonso.solve!(alf1) + @test Alfonso.get_niters(alf1) <= 55 + @test Alfonso.get_status(alf1) == :Optimal + + @test Alfonso.get_pobj(alf2) ≈ Alfonso.get_pobj(alf1) atol=1e-4 rtol=1e-4 + @test Alfonso.get_dobj(alf2) ≈ Alfonso.get_dobj(alf1) atol=1e-4 rtol=1e-4 end @testset "3D poly envelope example (sparse A)" begin diff --git a/test/runtests.jl b/test/runtests.jl index 1226c727e..90a769f1a 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2,65 +2,65 @@ using Alfonso using Test -# verbflag = true # Alfonso verbose option +verbflag = true # Alfonso verbose option -# # TODO interpolation tests -# -# # load optimizer builder functions from examples folder -# egs_dir = joinpath(@__DIR__, "../examples") -# include(joinpath(egs_dir, "envelope/envelope.jl")) -# include(joinpath(egs_dir, "lp/lp.jl")) -# include(joinpath(egs_dir, "namedpoly/namedpoly.jl")) -# -# # run native and MOI interfaces on examples -# include(joinpath(@__DIR__, "nativeexamples.jl")) - - -import MathOptInterface -const MOI = MathOptInterface -const MOIT = MOI.Test -const MOIB = MOI.Bridges -const MOIU = MOI.Utilities +# TODO interpolation tests -MOIU.@model(AlfonsoModelData, - (), - (MOI.EqualTo, MOI.GreaterThan, MOI.LessThan), - (MOI.Zeros, MOI.Nonnegatives, MOI.Nonpositives, MOI.SecondOrderCone, MOI.RotatedSecondOrderCone, MOI.PositiveSemidefiniteConeTriangle, MOI.ExponentialCone, MOI.PowerCone), - (), - (MOI.SingleVariable,), - (MOI.ScalarAffineFunction,), - (MOI.VectorOfVariables,), - (MOI.VectorAffineFunction,), - ) +# load optimizer builder functions from examples folder +egs_dir = joinpath(@__DIR__, "../examples") +include(joinpath(egs_dir, "envelope/envelope.jl")) +include(joinpath(egs_dir, "lp/lp.jl")) +include(joinpath(egs_dir, "namedpoly/namedpoly.jl")) -const optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Optimizer()) +# run native and MOI interfaces on examples +include(joinpath(@__DIR__, "nativeexamples.jl")) -const config = MOIT.TestConfig( - atol=1e-4, - rtol=1e-4, - solve=true, - query=true, - modify_lhs=true, - duals=false, - infeas_certificates=false, - ) - -# @testset "Continuous linear problems" begin -# MOIT.contlineartest(MOIB.SplitInterval{Float64}(optimizer), config) +# +# import MathOptInterface +# const MOI = MathOptInterface +# const MOIT = MOI.Test +# const MOIB = MOI.Bridges +# const MOIU = MOI.Utilities +# +# MOIU.@model(AlfonsoModelData, +# (), +# (MOI.EqualTo, MOI.GreaterThan, MOI.LessThan), +# (MOI.Zeros, MOI.Nonnegatives, MOI.Nonpositives, MOI.SecondOrderCone, MOI.RotatedSecondOrderCone, MOI.PositiveSemidefiniteConeTriangle, MOI.ExponentialCone, MOI.PowerCone), +# (), +# (MOI.SingleVariable,), +# (MOI.ScalarAffineFunction,), +# (MOI.VectorOfVariables,), +# (MOI.VectorAffineFunction,), +# ) +# +# const optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Optimizer()) +# +# const config = MOIT.TestConfig( +# atol=1e-4, +# rtol=1e-4, +# solve=true, +# query=true, +# modify_lhs=true, +# duals=false, +# infeas_certificates=false, +# ) +# +# # @testset "Continuous linear problems" begin +# # MOIT.contlineartest(MOIB.SplitInterval{Float64}(optimizer), config) +# # end +# +# @testset "Continuous conic problems" begin +# exclude = ["rootdet", "logdet"] +# MOIT.contconictest( +# MOIB.SquarePSD{Float64}( +# MOIB.GeoMean{Float64}( +# # MOIB.LogDet{Float64}( +# # MOIB.RootDet{Float64}( +# optimizer +# )), +# config, exclude) # end -@testset "Continuous conic problems" begin - exclude = ["rootdet", "logdet"] - MOIT.contconictest( - MOIB.SquarePSD{Float64}( - MOIB.GeoMean{Float64}( - # MOIB.LogDet{Float64}( - # MOIB.RootDet{Float64}( - optimizer - )), - config, exclude) -end - # # From 440617ba81afc2b258ebc0e82bce92039d176b63 Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Sat, 29 Sep 2018 15:20:32 -0400 Subject: [PATCH 06/10] prep to transition to new format in MOI functions --- src/Alfonso.jl | 6 +-- src/cone.jl | 19 +++++++ src/mathoptinterface.jl | 59 +++++++++++++-------- test/runtests.jl | 114 +++++++++++++++++++++------------------- 4 files changed, 120 insertions(+), 78 deletions(-) diff --git a/src/Alfonso.jl b/src/Alfonso.jl index 403b4f936..4ac770641 100644 --- a/src/Alfonso.jl +++ b/src/Alfonso.jl @@ -20,7 +20,7 @@ module Alfonso end include("nativeinterface.jl") - # import MathOptInterface - # const MOI = MathOptInterface - # include("mathoptinterface.jl") + import MathOptInterface + const MOI = MathOptInterface + include("mathoptinterface.jl") end diff --git a/src/cone.jl b/src/cone.jl index 8d941d266..295c09631 100644 --- a/src/cone.jl +++ b/src/cone.jl @@ -9,6 +9,25 @@ abstract type PrimitiveCone end mutable struct Cone prms::Vector{PrimitiveCone} idxs::Vector{AbstractVector{Int}} + + function Cone(prms::Vector{PrimitiveCone}, idxs::Vector{AbstractVector{Int}}) + @assert length(prms) == length(idxs) + for k in eachindex(prms) + @assert dimension(prms[k]) == length(idxs[k]) + end + cone = new() + cone.prms = prms + cone.idxs = idxs + return cone + end +end +Cone() = Cone(PrimitiveCone[], AbstractVector{Int}[]) + +function addprimitivecone!(cone::Cone, prm::PrimitiveCone, idx::AbstractVector{Int}) + @assert dimension(prm) == length(idx) + push!(cone.prms, prm) + push!(cone.idxs, idx) + return cone end # calculate complexity parameter of the barrier (sum of the primitive cone barrier parameters) diff --git a/src/mathoptinterface.jl b/src/mathoptinterface.jl index 8ba97f1d3..086d77a8b 100644 --- a/src/mathoptinterface.jl +++ b/src/mathoptinterface.jl @@ -44,6 +44,13 @@ const SupportedSets = Union{ MOI.supports_constraint(::Optimizer, ::Type{<:SupportedFuns}, ::Type{<:SupportedSets}) = true +conefrommoi(s::MOI.SecondOrderCone) = SecondOrderCone(MOI.dimension(s)) +conefrommoi(s::MOI.RotatedSecondOrderCone) = RotatedSecondOrderCone(MOI.dimension(s)) +conefrommoi(s::MOI.PositiveSemidefiniteConeTriangle) = PositiveSemidefiniteCone(MOI.dimension(s)) +conefrommoi(s::MOI.ExponentialCone) = ExponentialCone() +# conefrommoi(s::MOI.PowerCone) = PowerCone(s.exponent) + + # build representation as min c'x s.t. Ax = b, x in K # TODO what if some variables x are in multiple cone constraints? function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_attributes=true) @@ -79,15 +86,12 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ getconset(conidx) = MOI.get(src, MOI.ConstraintSet(), conidx) # pass over vector constraints to construct conic model data - (IA, JA, VA) = (Int[], Int[], Float64[]) - (Ib, Vb) = (Int[], Float64[]) - cone = Cone() - i = 0 # MOI constraint objects - p = 0 # rows of A (equality constraint matrix) - q = 0 # rows of G (cone constraint matrix) # equality constraints + p = 0 # rows of A (equality constraint matrix) + (IA, JA, VA) = (Int[], Int[], Float64[]) + (Ib, Vb) = (Int[], Float64[]) for ci in getsrccons(MOI.SingleVariable, MOI.EqualTo{Float64}) i += 1 @@ -169,17 +173,19 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ p += dim end + # conic constraints + q = 0 # rows of G (cone constraint matrix) + (IG, JG, VG) = (Int[], Int[], Float64[]) + (Ih, Vh) = (Int[], Float64[]) + cone = Cone() + # LP constraints: build up one nonnegative cone and/or one nonpositive cone + nonnegvars = Int[] + nonposvars = Int[] + # TODO also use variable bounds to build one L_inf cone - - - - - # linear inequality constraints - - nonnegvars = Int[] # for building up a single nonnegative cone - # TODO do the same for variable bounds to make a L_inf cone +# TODO use nonpositive cone for S in (MOI.GreaterThan{Float64}, MOI.LessThan{Float64}) for ci in getsrccons(MOI.SingleVariable, S) @@ -266,16 +272,19 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ end end + + + # add single nonnegative cone addprimitivecone!(cone, NonnegativeCone(length(nonnegvars)), nonnegvars) - # add other primitive cone constraints (in increasing order of difficulty of evaluating in-cone checks) - conefrommoi(s::MOI.SecondOrderCone) = SecondOrderCone(MOI.dimension(s)) - conefrommoi(s::MOI.RotatedSecondOrderCone) = RotatedSecondOrderCone(MOI.dimension(s)) - conefrommoi(s::MOI.PositiveSemidefiniteConeTriangle) = PositiveSemidefiniteCone(MOI.dimension(s)) - conefrommoi(s::MOI.ExponentialCone) = ExponentialCone() - # conefrommoi(s::MOI.PowerCone) = PowerCone(s.exponent) + + + + + + # add other primitive cone constraints (in increasing order of difficulty of evaluating in-cone checks) for S in (MOI.SecondOrderCone, MOI.RotatedSecondOrderCone, MOI.ExponentialCone, MOI.PowerCone, MOI.PositiveSemidefiniteConeTriangle) for ci in getsrccons(MOI.VectorOfVariables, S) @@ -314,11 +323,17 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ end end + + + + # finalize data and load into Alfonso model + c = Vector(dropzeros!(sparsevec(Jc, Vc, n))) A = dropzeros!(sparse(IA, JA, VA, p, n)) b = Vector(dropzeros!(sparsevec(Ib, Vb, p))) - c = Vector(dropzeros!(sparsevec(Jc, Vc, n))) - load_data!(opt.alf, A, b, c, cone) + G = dropzeros!(sparse(IG, JG, VG, q, n)) + h = Vector(dropzeros!(sparsevec(Ih, Vh, q))) + Alfonso.load_data!(alf, c, A, b, G, h, cone) # set objective sense opt.objsense = MOI.get(src, MOI.ObjectiveSense()) diff --git a/test/runtests.jl b/test/runtests.jl index 90a769f1a..6b5ddaf56 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2,65 +2,73 @@ using Alfonso using Test -verbflag = true # Alfonso verbose option +# verbflag = true # Alfonso verbose option +# +# # TODO interpolation tests +# +# # load optimizer builder functions from examples folder +# egs_dir = joinpath(@__DIR__, "../examples") +# include(joinpath(egs_dir, "envelope/envelope.jl")) +# include(joinpath(egs_dir, "lp/lp.jl")) +# include(joinpath(egs_dir, "namedpoly/namedpoly.jl")) +# +# # run native and MOI interfaces on examples +# include(joinpath(@__DIR__, "nativeexamples.jl")) -# TODO interpolation tests -# load optimizer builder functions from examples folder -egs_dir = joinpath(@__DIR__, "../examples") -include(joinpath(egs_dir, "envelope/envelope.jl")) -include(joinpath(egs_dir, "lp/lp.jl")) -include(joinpath(egs_dir, "namedpoly/namedpoly.jl")) +import MathOptInterface +const MOI = MathOptInterface +const MOIT = MOI.Test +const MOIB = MOI.Bridges +const MOIU = MOI.Utilities -# run native and MOI interfaces on examples -include(joinpath(@__DIR__, "nativeexamples.jl")) +MOIU.@model(AlfonsoModelData, + (), + ( + MOI.EqualTo, MOI.GreaterThan, MOI.LessThan, + # MOI.Interval, + ), + ( + MOI.Zeros, MOI.Nonnegatives, MOI.Nonpositives, + MOI.SecondOrderCone, MOI.RotatedSecondOrderCone, + MOI.PositiveSemidefiniteConeTriangle, + MOI.ExponentialCone, + # MOI.PowerCone, + ), + (), + (MOI.SingleVariable,), + (MOI.ScalarAffineFunction,), + (MOI.VectorOfVariables,), + (MOI.VectorAffineFunction,), + ) -# -# import MathOptInterface -# const MOI = MathOptInterface -# const MOIT = MOI.Test -# const MOIB = MOI.Bridges -# const MOIU = MOI.Utilities -# -# MOIU.@model(AlfonsoModelData, -# (), -# (MOI.EqualTo, MOI.GreaterThan, MOI.LessThan), -# (MOI.Zeros, MOI.Nonnegatives, MOI.Nonpositives, MOI.SecondOrderCone, MOI.RotatedSecondOrderCone, MOI.PositiveSemidefiniteConeTriangle, MOI.ExponentialCone, MOI.PowerCone), -# (), -# (MOI.SingleVariable,), -# (MOI.ScalarAffineFunction,), -# (MOI.VectorOfVariables,), -# (MOI.VectorAffineFunction,), -# ) -# -# const optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Optimizer()) -# -# const config = MOIT.TestConfig( -# atol=1e-4, -# rtol=1e-4, -# solve=true, -# query=true, -# modify_lhs=true, -# duals=false, -# infeas_certificates=false, -# ) -# -# # @testset "Continuous linear problems" begin -# # MOIT.contlineartest(MOIB.SplitInterval{Float64}(optimizer), config) -# # end -# -# @testset "Continuous conic problems" begin -# exclude = ["rootdet", "logdet"] -# MOIT.contconictest( -# MOIB.SquarePSD{Float64}( -# MOIB.GeoMean{Float64}( -# # MOIB.LogDet{Float64}( -# # MOIB.RootDet{Float64}( -# optimizer -# )), -# config, exclude) +const optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Optimizer()) + +const config = MOIT.TestConfig( + atol=1e-4, + rtol=1e-4, + solve=true, + query=true, + modify_lhs=true, + duals=false, + infeas_certificates=false, + ) + +# @testset "Continuous linear problems" begin +# MOIT.contlineartest(MOIB.SplitInterval{Float64}(optimizer), config) # end +@testset "Continuous conic problems" begin + MOIT.contconictest( + MOIB.SquarePSD{Float64}( + MOIB.GeoMean{Float64}( + MOIB.LogDet{Float64}( + MOIB.RootDet{Float64}( + optimizer + )))), + config) +end + # # From 902371a5c8dd59db64971d2f1e8cbfee78943b41 Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Sat, 29 Sep 2018 16:51:34 -0400 Subject: [PATCH 07/10] construct model in copy_to --- Manifest.toml | 6 +- src/Alfonso.jl | 2 +- src/mathoptinterface.jl | 150 +++++++++++++++++----------------------- test/runtests.jl | 17 ++--- 4 files changed, 77 insertions(+), 98 deletions(-) diff --git a/Manifest.toml b/Manifest.toml index 5460bf30c..e066ecad1 100644 --- a/Manifest.toml +++ b/Manifest.toml @@ -79,8 +79,10 @@ uuid = "d6f4376e-aef5-505a-96c1-9c027394607a" [[MathOptInterface]] deps = ["Compat", "Unicode"] git-tree-sha1 = "ba12e7ce825c1458c03f88aae84fa630d882d303" +repo-rev = "master" +repo-url = "https://github.com/JuliaOpt/MathOptInterface.jl.git" uuid = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" -version = "0.6.1" +version = "0.6.1+" [[Mmap]] uuid = "a63ad114-7e13-5084-954f-fe012c677804" @@ -139,7 +141,7 @@ deps = ["Distributed", "InteractiveUtils", "Logging", "Random"] uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [[UUIDs]] -deps = ["Random", "SHA"] +deps = ["Random"] uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4" [[Unicode]] diff --git a/src/Alfonso.jl b/src/Alfonso.jl index 4ac770641..2c00bd4f6 100644 --- a/src/Alfonso.jl +++ b/src/Alfonso.jl @@ -21,6 +21,6 @@ module Alfonso include("nativeinterface.jl") import MathOptInterface - const MOI = MathOptInterface + MOI = MathOptInterface include("mathoptinterface.jl") end diff --git a/src/mathoptinterface.jl b/src/mathoptinterface.jl index 086d77a8b..5f6a19dcb 100644 --- a/src/mathoptinterface.jl +++ b/src/mathoptinterface.jl @@ -163,7 +163,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ fi = getconfun(ci) dim = MOI.output_dimension(fi) - append!(Ib, collect(p+1:p+dim)) + append!(Ib, p+1:p+dim) append!(Vb, -fi.constants) for vt in fi.terms push!(IA, p + vt.output_index) @@ -180,36 +180,27 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ cone = Cone() # LP constraints: build up one nonnegative cone and/or one nonpositive cone - nonnegvars = Int[] - nonposvars = Int[] - # TODO also use variable bounds to build one L_inf cone - + nonnegrows = Int[] + nonposrows = Int[] -# TODO use nonpositive cone + # TODO also use variable bounds to build one L_inf cone for S in (MOI.GreaterThan{Float64}, MOI.LessThan{Float64}) for ci in getsrccons(MOI.SingleVariable, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, S}(i) - bval = ((S == MOI.GreaterThan{Float64}) ? getconset(ci).lower : getconset(ci).upper) - vn = getconfun(ci).variable - if iszero(bval) - # no need to change A,b - push!(nonnegvars, idxmap[vn].value) + q += 1 + push!(IG, q) + push!(JG, idxmap[getconfun(ci).variable].value) + push!(VG, -1.0) + push!(Ih, q) + if S == MOI.GreaterThan{Float64} + push!(Vh, -getconset(ci).lower) + push!(nonnegrows, q) else - # add auxiliary variable and equality constraint - p += 1 - push!(Ib, p) - push!(Vb, bval) - push!(IA, p) - push!(JA, idxmap[vn].value) - push!(VA, 1.0) - n += 1 - push!(IA, p) - push!(JA, n) - push!(VA, ((S == MOI.GreaterThan{Float64}) ? -1.0 : 1.0)) - push!(nonnegvars, n) + push!(Vh, -getconset(ci).upper) + push!(nonposrows, q) end end @@ -217,21 +208,21 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, S}(i) - # add auxiliary variable and equality constraint fi = getconfun(ci) - p += 1 - push!(Ib, p) - push!(Vb, ((S == MOI.GreaterThan{Float64}) ? getconset(ci).lower : getconset(ci).upper) - fi.constant) + q += 1 for vt in fi.terms - push!(IA, p) - push!(JA, idxmap[vt.variable_index].value) - push!(VA, vt.coefficient) + push!(IG, q) + push!(JG, idxmap[vt.variable_index].value) + push!(VG, -vt.coefficient) + end + push!(Ih, q) + if S == MOI.GreaterThan{Float64} + push!(Vh, -getconset(ci).lower + fi.constant) + push!(nonnegrows, q) + else + push!(Vh, -getconset(ci).upper + fi.constant) + push!(nonposrows, q) end - n += 1 - push!(IA, p) - push!(JA, n) - push!(VA, ((S == MOI.GreaterThan{Float64}) ? -1.0 : 1.0)) - push!(nonnegvars, n) end end @@ -241,12 +232,15 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, S}(i) for vj in getconfun(ci).variables - col = idxmap[vj].value - p += 1 - push!(IA, p) - push!(JA, col) - push!(VA, ((S == MOI.Nonnegatives) ? 1.0 : -1.0)) - push!(nonnegvars, col) + q += 1 + push!(IG, q) + push!(JG, idxmap[vj].value) + push!(VG, -1.0) + if S == MOI.Nonnegatives + push!(nonnegrows, q) + else + push!(nonposrows, q) + end end end @@ -257,34 +251,25 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ fi = getconfun(ci) dim = MOI.output_dimension(fi) for vt in fi.terms - push!(IA, p + vt.output_index) - push!(JA, idxmap[vt.scalar_term.variable_index].value) - push!(VA, vt.scalar_term.coefficient) + push!(IG, q + vt.output_index) + push!(JG, idxmap[vt.scalar_term.variable_index].value) + push!(VG, -vt.scalar_term.coefficient) + end + append!(Ih, q+1:q+dim) + append!(Vh, fi.constants) + if S == MOI.Nonnegatives + append!(nonnegrows, q+1:q+dim) + else + append!(nonposrows, q+1:q+dim) end - append!(Ib, collect(p+1:p+dim)) - append!(Vb, -fi.constants) - append!(IA, collect(p+1:p+dim)) - append!(JA, collect(n+1:n+dim)) - append!(VA, fill(((S == MOI.Nonnegatives) ? 1.0 : -1.0), dim)) - p += dim - append!(nonnegvars, collect(n+1:n+dim)) - n += dim + q += dim end end + addprimitivecone!(cone, NonnegativeCone(length(nonnegrows)), nonnegrows) + addprimitivecone!(cone, NonpositiveCone(length(nonposrows)), nonposrows) - - - # add single nonnegative cone - addprimitivecone!(cone, NonnegativeCone(length(nonnegvars)), nonnegvars) - - - - - - - - # add other primitive cone constraints (in increasing order of difficulty of evaluating in-cone checks) + # non-LP conic constraints for S in (MOI.SecondOrderCone, MOI.RotatedSecondOrderCone, MOI.ExponentialCone, MOI.PowerCone, MOI.PositiveSemidefiniteConeTriangle) for ci in getsrccons(MOI.VectorOfVariables, S) @@ -292,13 +277,12 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, S}(i) fi = getconfun(ci) - for vj in fi.variables - p += 1 - push!(IA, p) - push!(JA, idxmap[vj].value) - push!(VA, 1.0) - end - addprimitivecone!(cone, conefrommoi(getconset(ci)), [idxmap[vj].value for vj in fi.variables]) + dim = MOI.output_dimension(fi) + append!(IG, q+1:q+dim) + append!(JG, idxmap[vj].value for vj in fi.variables) + append!(VG, -ones(dim)) + addprimitivecone!(cone, conefrommoi(getconset(ci)), q+1:q+dim) + q += dim end for ci in getsrccons(MOI.VectorAffineFunction{Float64}, S) @@ -308,32 +292,24 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ fi = getconfun(ci) dim = MOI.output_dimension(fi) for vt in fi.terms - push!(IA, p + vt.output_index) - push!(JA, idxmap[vt.scalar_term.variable_index].value) - push!(VA, vt.scalar_term.coefficient) + push!(IG, q + vt.output_index) + push!(JG, idxmap[vt.scalar_term.variable_index].value) + push!(VG, -vt.scalar_term.coefficient) end - append!(Ib, collect(p+1:p+dim)) - append!(Vb, -fi.constants) - append!(IA, collect(p+1:p+dim)) - append!(JA, collect(n+1:n+dim)) - append!(VA, fill(1.0, dim)) - p += dim - addprimitivecone!(cone, conefrommoi(getconset(ci)), n+1:n+dim) - n += dim + append!(Ih, q+1:q+dim) + append!(Vh, fi.constants) + addprimitivecone!(cone, conefrommoi(getconset(ci)), q+1:q+dim) + q += dim end end - - - - # finalize data and load into Alfonso model c = Vector(dropzeros!(sparsevec(Jc, Vc, n))) A = dropzeros!(sparse(IA, JA, VA, p, n)) b = Vector(dropzeros!(sparsevec(Ib, Vb, p))) G = dropzeros!(sparse(IG, JG, VG, q, n)) h = Vector(dropzeros!(sparsevec(Ih, Vh, q))) - Alfonso.load_data!(alf, c, A, b, G, h, cone) + Alfonso.load_data!(opt.alf, c, A, b, G, h, cone) # set objective sense opt.objsense = MOI.get(src, MOI.ObjectiveSense()) diff --git a/test/runtests.jl b/test/runtests.jl index 6b5ddaf56..ed5441446 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -17,10 +17,10 @@ using Test import MathOptInterface -const MOI = MathOptInterface -const MOIT = MOI.Test -const MOIB = MOI.Bridges -const MOIU = MOI.Utilities +MOI = MathOptInterface +MOIT = MOI.Test +MOIB = MOI.Bridges +MOIU = MOI.Utilities MOIU.@model(AlfonsoModelData, (), @@ -42,13 +42,13 @@ MOIU.@model(AlfonsoModelData, (MOI.VectorAffineFunction,), ) -const optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Optimizer()) +optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Optimizer()) -const config = MOIT.TestConfig( +config = MOIT.TestConfig( atol=1e-4, rtol=1e-4, solve=true, - query=true, + query=false, modify_lhs=true, duals=false, infeas_certificates=false, @@ -59,6 +59,7 @@ const config = MOIT.TestConfig( # end @testset "Continuous conic problems" begin + exclude = ["rootdet", "logdet"] # TODO bridges not working? should not need to exclude in future MOIT.contconictest( MOIB.SquarePSD{Float64}( MOIB.GeoMean{Float64}( @@ -66,7 +67,7 @@ const config = MOIT.TestConfig( MOIB.RootDet{Float64}( optimizer )))), - config) + config, exclude) end From 1249639f51306d702a35d538cad76d11acd125f7 Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Sat, 29 Sep 2018 19:22:10 -0400 Subject: [PATCH 08/10] progress in making MOI getter functions work --- src/mathoptinterface.jl | 288 ++++++++++++++++++---------------------- test/runtests.jl | 199 ++++++++++++++------------- 2 files changed, 225 insertions(+), 262 deletions(-) diff --git a/src/mathoptinterface.jl b/src/mathoptinterface.jl index 5f6a19dcb..342f79d9b 100644 --- a/src/mathoptinterface.jl +++ b/src/mathoptinterface.jl @@ -1,7 +1,19 @@ mutable struct Optimizer <: MOI.AbstractOptimizer alf::AlfonsoOpt + objsense::MOI.OptimizationSense + constroffsets::Vector{Int} + numeqconstrs::Int + + x + s + y + z + status + solvetime + pobj + dobj function Optimizer(alf::AlfonsoOpt) opt = new() @@ -55,7 +67,7 @@ conefrommoi(s::MOI.ExponentialCone) = ExponentialCone() # TODO what if some variables x are in multiple cone constraints? function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_attributes=true) @assert !copy_names - idxmap = Dict{MOI.Index,MOI.Index}() + idxmap = Dict{MOI.Index, MOI.Index}() # model # mattr_src = MOI.get(src, MOI.ListOfModelAttributesSet()) @@ -82,10 +94,11 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ end getsrccons(F, S) = MOI.get(src, MOI.ListOfConstraintIndices{F,S}()) - getconfun(conidx) = MOI.get(src, MOI.ConstraintFunction(), conidx) - getconset(conidx) = MOI.get(src, MOI.ConstraintSet(), conidx) + getconfun(conidx) = MOI.get(src, MOI.ConstraintFunction(), conidx) + getconset(conidx) = MOI.get(src, MOI.ConstraintSet(), conidx) # pass over vector constraints to construct conic model data + constroffsets = Vector{Int}() i = 0 # MOI constraint objects # equality constraints @@ -93,86 +106,76 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ (IA, JA, VA) = (Int[], Int[], Float64[]) (Ib, Vb) = (Int[], Float64[]) - for ci in getsrccons(MOI.SingleVariable, MOI.EqualTo{Float64}) - i += 1 - idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, MOI.EqualTo{Float64}}(i) - - p += 1 - push!(Ib, p) - push!(Vb, getconset(ci).value) - push!(IA, p) - push!(JA, idxmap[getconfun(ci).variable].value) - push!(VA, 1.0) - end - - for ci in getsrccons(MOI.SingleVariable, MOI.Zeros) - i += 1 - idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, MOI.Zeros}(i) - - p += 1 - push!(IA, p) - push!(JA, idxmap[getconfun(ci).variable].value) - push!(VA, 1.0) - end - - for ci in getsrccons(MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}) - i += 1 - idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}}(i) + for S in (MOI.EqualTo{Float64}, MOI.Zeros) + # TODO can preprocess variables equal to constant + for ci in getsrccons(MOI.SingleVariable, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, S}(i) + push!(constroffsets, p) - fi = getconfun(ci) - p += 1 - push!(Ib, p) - push!(Vb, getconset(ci).value - fi.constant) - for vt in fi.terms + p += 1 push!(IA, p) - push!(JA, idxmap[vt.variable_index].value) - push!(VA, vt.coefficient) + push!(JA, idxmap[getconfun(ci).variable].value) + push!(VA, 1.0) + if S == MOI.EqualTo{Float64} + push!(Ib, p) + push!(Vb, getconset(ci).value) + end end - end - for ci in getsrccons(MOI.ScalarAffineFunction{Float64}, MOI.Zeros) - i += 1 - idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, MOI.Zeros}(i) + for ci in getsrccons(MOI.ScalarAffineFunction{Float64}, S) + i += 1 + idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, S}(i) + push!(constroffsets, p) - fi = getconfun(ci) - p += 1 - push!(Ib, p) - push!(Vb, -fi.constant) - for vt in fi.terms - push!(IA, p) - push!(JA, idxmap[vt.variable_index].value) - push!(VA, vt.coefficient) + p += 1 + fi = getconfun(ci) + for vt in fi.terms + push!(IA, p) + push!(JA, idxmap[vt.variable_index].value) + push!(VA, vt.coefficient) + end + push!(Ib, p) + if S == MOI.EqualTo{Float64} + push!(Vb, getconset(ci).value - fi.constant) + else + push!(Vb, -fi.constant) + end end end + # TODO can preprocess variables equal to zero for ci in getsrccons(MOI.VectorOfVariables, MOI.Zeros) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, MOI.Zeros}(i) + push!(constroffsets, p) - for vi in getconfun(ci).variables - p += 1 - push!(IA, p) - push!(JA, idxmap[vi].value) - push!(VA, 1.0) - end + fi = getconfun(ci) + dim = MOI.output_dimension(fi) + append!(IA, p+1:p+dim) + append!(JA, idxmap[vi].value for vi in fi.variables) + append!(VA, ones(dim)) end for ci in getsrccons(MOI.VectorAffineFunction{Float64}, MOI.Zeros) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, MOI.Zeros}(i) + push!(constroffsets, p) fi = getconfun(ci) dim = MOI.output_dimension(fi) - append!(Ib, p+1:p+dim) - append!(Vb, -fi.constants) for vt in fi.terms push!(IA, p + vt.output_index) push!(JA, idxmap[vt.scalar_term.variable_index].value) push!(VA, vt.scalar_term.coefficient) end + append!(Ib, p+1:p+dim) + append!(Vb, -fi.constants) p += dim end + numeqconstrs = i + # conic constraints q = 0 # rows of G (cone constraint matrix) (IG, JG, VG) = (Int[], Int[], Float64[]) @@ -189,6 +192,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.SingleVariable, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, S}(i) + push!(constroffsets, q) q += 1 push!(IG, q) @@ -207,9 +211,10 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.ScalarAffineFunction{Float64}, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, S}(i) + push!(constroffsets, q) - fi = getconfun(ci) q += 1 + fi = getconfun(ci) for vt in fi.terms push!(IG, q) push!(JG, idxmap[vt.variable_index].value) @@ -230,6 +235,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.VectorOfVariables, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, S}(i) + push!(constroffsets, q) for vj in getconfun(ci).variables q += 1 @@ -247,6 +253,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.VectorAffineFunction{Float64}, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, S}(i) + push!(constroffsets, q) fi = getconfun(ci) dim = MOI.output_dimension(fi) @@ -275,6 +282,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.VectorOfVariables, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, S}(i) + push!(constroffsets, q) fi = getconfun(ci) dim = MOI.output_dimension(fi) @@ -288,6 +296,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.VectorAffineFunction{Float64}, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, S}(i) + push!(constroffsets, q) fi = getconfun(ci) dim = MOI.output_dimension(fi) @@ -309,26 +318,42 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ b = Vector(dropzeros!(sparsevec(Ib, Vb, p))) G = dropzeros!(sparse(IG, JG, VG, q, n)) h = Vector(dropzeros!(sparsevec(Ih, Vh, q))) - Alfonso.load_data!(opt.alf, c, A, b, G, h, cone) + # TODO make converting to dense optional + # Alfonso.load_data!(opt.alf, c, A, b, G, h, cone) # sparse + Alfonso.load_data!(opt.alf, c, Matrix(A), b, Matrix(G), h, cone) # dense - # set objective sense + # store information needed by MOI getter functions opt.objsense = MOI.get(src, MOI.ObjectiveSense()) + opt.constroffsets = constroffsets + opt.numeqconstrs = numeqconstrs return idxmap end -MOI.optimize!(opt::Optimizer) = solve!(opt.alf) +function MOI.optimize!(opt::Optimizer) + solve!(opt.alf) + + opt.x = get_x(opt.alf) + opt.s = get_s(opt.alf) + opt.y = get_y(opt.alf) + opt.z = get_z(opt.alf) + opt.status = get_status(opt.alf) + opt.solvetime = get_solvetime(opt.alf) + opt.pobj = get_pobj(opt.alf) + opt.dobj = get_dobj(opt.alf) + + return nothing +end # function MOI.free!(opt::Optimizer) # TODO call gc on opt.alf? function MOI.get(opt::Optimizer, ::MOI.TerminationStatus) # TODO time limit etc - status = get_status(opt.alf) - if status in (:Optimal, :NearlyInfeasible, :IllPosed) + if opt.status in (:Optimal, :PrimalInfeasible, :DualInfeasible, :IllPosed) return MOI.Success - elseif status in (:PredictorFail, :CorrectorFail) + elseif opt.status in (:PredictorFail, :CorrectorFail) return MOI.NumericalError - elseif status == :IterationLimit + elseif opt.status == :IterationLimit return MOI.IterationLimit else return MOI.OtherError @@ -337,116 +362,79 @@ end function MOI.get(opt::Optimizer, ::MOI.ObjectiveValue) if opt.objsense == MOI.MinSense - return get_pobj(opt.alf) + return opt.pobj elseif opt.objsense == MOI.MaxSense - return -get_pobj(opt.alf) + return -opt.pobj else - return NaN + error("no objective sense is set") end end function MOI.get(opt::Optimizer, ::MOI.ObjectiveBound) if opt.objsense == MOI.MinSense - return get_dobj(opt.alf) + return opt.dobj elseif opt.objsense == MOI.MaxSense - return -get_dobj(opt.alf) + return -opt.dobj else - return NaN + error("no objective sense is set") end end function MOI.get(opt::Optimizer, ::MOI.ResultCount) - status = get_status(opt.alf) - if status in (:Optimal, :NearlyInfeasible, :IllPosed) + if opt.status in (:Optimal, :PrimalInfeasible, :DualInfeasible, :IllPosed) return 1 end return 0 end function MOI.get(opt::Optimizer, ::MOI.PrimalStatus) - status = get_status(opt.alf) - if status == :Optimal + if opt.status == :Optimal return MOI.FeasiblePoint + elseif opt.status == :PrimalInfeasible + return MOI.InfeasiblePoint + elseif opt.status == :DualInfeasible + return MOI.InfeasibilityCertificate + elseif opt.status == :IllPosed + return MOI.UnknownResultStatus # TODO later distinguish primal/dual ill posed certificates else return MOI.UnknownResultStatus end end function MOI.get(opt::Optimizer, ::MOI.DualStatus) - status = get_status(opt.alf) - if status == :Optimal + if opt.status == :Optimal return MOI.FeasiblePoint + elseif opt.status == :PrimalInfeasible + return MOI.InfeasibilityCertificate + elseif opt.status == :DualInfeasible + return MOI.InfeasiblePoint + elseif opt.status == :IllPosed + return MOI.UnknownResultStatus # TODO later distinguish primal/dual ill posed certificates else return MOI.UnknownResultStatus end end +MOI.get(opt::Optimizer, ::MOI.VariablePrimal, vi::MOI.VariableIndex) = opt.x[vi.value] +MOI.get(opt::Optimizer, a::MOI.VariablePrimal, vi::Vector{MOI.VariableIndex}) = MOI.get.(opt, a, vi) -# TODO from ECOS: -# Implements getter for result value and statuses -# function MOI.get(instance::Optimizer, ::MOI.TerminationStatus) -# flag = instance.sol.ret_val -# if flag == ECOS.ECOS_OPTIMAL -# MOI.Success -# elseif flag == ECOS.ECOS_PINF -# MOI.Success -# elseif flag == ECOS.ECOS_DINF # Dual infeasible = primal unbounded, probably -# MOI.Success -# elseif flag == ECOS.ECOS_MAXIT -# MOI.IterationLimit -# elseif flag == ECOS.ECOS_OPTIMAL + ECOS.ECOS_INACC_OFFSET -# m.solve_stat = MOI.AlmostSuccess -# else -# m.solve_stat = MOI.OtherError -# end -# end -# -# MOI.get(instance::Optimizer, ::MOI.ObjectiveValue) = instance.sol.objval -# MOI.get(instance::Optimizer, ::MOI.ObjectiveBound) = instance.sol.objbnd -# -# function MOI.get(instance::Optimizer, ::MOI.PrimalStatus) -# flag = instance.sol.ret_val -# if flag == ECOS.ECOS_OPTIMAL -# MOI.FeasiblePoint -# elseif flag == ECOS.ECOS_PINF -# MOI.InfeasiblePoint -# elseif flag == ECOS.ECOS_DINF # Dual infeasible = primal unbounded, probably -# MOI.InfeasibilityCertificate -# elseif flag == ECOS.ECOS_MAXIT -# MOI.UnknownResultStatus -# elseif flag == ECOS.ECOS_OPTIMAL + ECOS.ECOS_INACC_OFFSET -# m.solve_stat = MOI.NearlyFeasiblePoint -# else -# m.solve_stat = MOI.OtherResultStatus -# end -# end -# function MOI.get(instance::Optimizer, ::MOI.DualStatus) -# flag = instance.sol.ret_val -# if flag == ECOS.ECOS_OPTIMAL -# MOI.FeasiblePoint -# elseif flag == ECOS.ECOS_PINF -# MOI.InfeasibilityCertificate -# elseif flag == ECOS.ECOS_DINF # Dual infeasible = primal unbounded, probably -# MOI.InfeasiblePoint -# elseif flag == ECOS.ECOS_MAXIT -# MOI.UnknownResultStatus -# elseif flag == ECOS.ECOS_OPTIMAL + ECOS.ECOS_INACC_OFFSET -# m.solve_stat = MOI.NearlyFeasiblePoint -# else -# m.solve_stat = MOI.OtherResultStatus -# end -# end - - -MOI.get(opt::Optimizer, ::MOI.VariablePrimal, vi::MOI.VariableIndex) = get_x(opt.alf)[vi.value] -MOI.get(opt::Optimizer, a::MOI.VariablePrimal, vi::Vector{MOI.VariableIndex}) = MOI.get.(opt, Ref(a), vi) +function MOI.get(opt::Optimizer, ::MOI.ConstraintDual, ci::MOI.ConstraintIndex{F, S}) where {F, S} + os = opt.constroffsets + i = ci.value + if i <= opt.numeqconstrs + # constraint is an equality + return opt.y[os[i]+1:os[i+1]] + else + # constraint is conic + return opt.z[os[i]+1:os[i+1]] + end +end +MOI.get(opt::Optimizer, a::MOI.ConstraintDual, ci::Vector{MOI.ConstraintIndex}) = MOI.get.(opt, a, ci) -# function MOI.get(opt::Optimizer, ::MOI.ConstraintPrimal, ci::MOI.ConstraintIndex{<:MOI.AbstractFunction, MOI.AbstractSet}) -# end -# MOI.get(opt::Optimizer, a::MOI.ConstraintPrimal, ci::Vector{MOI.ConstraintIndex}) = MOI.get.(opt, Ref(a), ci) +# TODO use x and s. if conic then just return s. otherwise have to use x - rhs after saving rhs in opt maybe +MOI.get(opt::Optimizer, ::MOI.ConstraintPrimal, ci::MOI.ConstraintIndex) = NaN +MOI.get(opt::Optimizer, a::MOI.ConstraintPrimal, ci::Vector{MOI.ConstraintIndex}) = MOI.get.(opt, a, ci) -MOI.get(opt::Optimizer, ::MOI.ConstraintDual, ci::MOI.ConstraintIndex) = get_y(opt.alf)[ci.value] -MOI.get(opt::Optimizer, a::MOI.ConstraintDual, ci::Vector{MOI.ConstraintIndex}) = MOI.get.(opt, Ref(a), ci) @@ -474,27 +462,3 @@ MOI.get(opt::Optimizer, a::MOI.ConstraintDual, ci::Vector{MOI.ConstraintIndex}) # end # # MOI.get(instance::Optimizer, ::MOI.ResultCount) = 1 - - - -# -# get_status(alf::AlfonsoOpt) = alf.status -# get_solvetime(alf::AlfonsoOpt) = alf.solvetime -# get_niters(alf::AlfonsoOpt) = alf.niters -# get_y(alf::AlfonsoOpt) = copy(alf.y) -# get_x(alf::AlfonsoOpt) = copy(alf.x) -# get_tau(alf::AlfonsoOpt) = alf.tau -# get_s(alf::AlfonsoOpt) = copy(alf.s) -# get_kappa(alf::AlfonsoOpt) = alf.kappa -# get_pobj(alf::AlfonsoOpt) = alf.pobj -# get_dobj(alf::AlfonsoOpt) = alf.dobj -# get_dgap(alf::AlfonsoOpt) = alf.dgap -# get_cgap(alf::AlfonsoOpt) = alf.cgap -# get_rel_dgap(alf::AlfonsoOpt) = alf.rel_dgap -# get_rel_cgap(alf::AlfonsoOpt) = alf.rel_cgap -# get_pres(alf::AlfonsoOpt) = copy(alf.pres) -# get_dres(alf::AlfonsoOpt) = copy(alf.dres) -# get_pin(alf::AlfonsoOpt) = alf.pin -# get_din(alf::AlfonsoOpt) = alf.din -# get_rel_pin(alf::AlfonsoOpt) = alf.rel_pin -# get_rel_din(alf::AlfonsoOpt) = alf.rel_din diff --git a/test/runtests.jl b/test/runtests.jl index ed5441446..7d42971c5 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -48,115 +48,114 @@ config = MOIT.TestConfig( atol=1e-4, rtol=1e-4, solve=true, - query=false, + query=true, modify_lhs=true, - duals=false, - infeas_certificates=false, + duals=true, + infeas_certificates=true, ) # @testset "Continuous linear problems" begin -# MOIT.contlineartest(MOIB.SplitInterval{Float64}(optimizer), config) +# MOIT.contlineartest( +# MOIB.SplitInterval{Float64}( +# optimizer +# ), +# config) # end -@testset "Continuous conic problems" begin - exclude = ["rootdet", "logdet"] # TODO bridges not working? should not need to exclude in future - MOIT.contconictest( - MOIB.SquarePSD{Float64}( - MOIB.GeoMean{Float64}( - MOIB.LogDet{Float64}( - MOIB.RootDet{Float64}( - optimizer - )))), - config, exclude) +# @testset "Continuous conic problems" begin +# exclude = ["rootdet", "logdet"] # TODO bridges not working? should not need to exclude in future +# MOIT.contconictest( +# MOIB.SquarePSD{Float64}( +# MOIB.GeoMean{Float64}( +# # MOIB.LogDet{Float64}( +# # MOIB.RootDet{Float64}( +# optimizer +# )),#)), +# config, exclude) +# end + + + +model = optimizer + +atol = config.atol +rtol = config.rtol + +# simple 2 variable, 1 constraint problem +# min -x +# st x + y <= 1 (x + y - 1 ∈ Nonpositives) +# x, y >= 0 (x, y ∈ Nonnegatives) + +@test MOI.supports(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) +@test MOI.supports(model, MOI.ObjectiveSense()) +@test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}) +@test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.GreaterThan{Float64}) +@test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}) +@test MOI.supports_constraint(model, MOI.SingleVariable, MOI.EqualTo{Float64}) +@test MOI.supports_constraint(model, MOI.SingleVariable, MOI.GreaterThan{Float64}) + +MOI.empty!(model) +@test MOI.is_empty(model) + +v = MOI.add_variables(model, 2) +@test MOI.get(model, MOI.NumberOfVariables()) == 2 + +cf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0,1.0], v), 0.0) +c = MOI.add_constraint(model, cf, MOI.LessThan(1.0)) +@test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}}()) == 1 + +vc1 = MOI.add_constraint(model, MOI.SingleVariable(v[1]), MOI.GreaterThan(0.0)) +# test fallback +vc2 = MOI.add_constraint(model, v[2], MOI.GreaterThan(0.0)) +@test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 2 + +# note: adding some redundant zero coefficients to catch solvers that don't handle duplicate coefficients correctly: +objf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([0.0,0.0,-1.0,0.0,0.0,0.0], [v; v; v]), 0.0) +MOI.set(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), objf) +MOI.set(model, MOI.ObjectiveSense(), MOI.MinSense) + +@test MOI.get(model, MOI.ObjectiveSense()) == MOI.MinSense + +if config.query + vrs = MOI.get(model, MOI.ListOfVariableIndices()) + @test vrs == v || vrs == reverse(v) + + @test objf ≈ MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) + + @test cf ≈ MOI.get(model, MOI.ConstraintFunction(), c) + + s = MOI.get(model, MOI.ConstraintSet(), c) + @test s == MOI.LessThan(1.0) + + s = MOI.get(model, MOI.ConstraintSet(), vc1) + @test s == MOI.GreaterThan(0.0) + + s = MOI.get(model, MOI.ConstraintSet(), vc2) + @test s == MOI.GreaterThan(0.0) end +if config.solve + MOI.optimize!(model) -# -# -# -# model = optimizer -# -# atol = config.atol -# rtol = config.rtol -# # simple 2 variable, 1 constraint problem -# # min -x -# # st x + y <= 1 (x + y - 1 ∈ Nonpositives) -# # x, y >= 0 (x, y ∈ Nonnegatives) -# -# @test MOI.supports(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) -# @test MOI.supports(model, MOI.ObjectiveSense()) -# @test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}) -# @test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.GreaterThan{Float64}) -# @test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}) -# @test MOI.supports_constraint(model, MOI.SingleVariable, MOI.EqualTo{Float64}) -# @test MOI.supports_constraint(model, MOI.SingleVariable, MOI.GreaterThan{Float64}) -# -# #@test MOI.get(model, MOI.SupportsAddConstraintAfterSolve()) -# #@test MOI.get(model, MOI.SupportsAddVariableAfterSolve()) -# #@test MOI.get(model, MOI.SupportsDeleteConstraint()) -# -# MOI.empty!(model) -# @test MOI.is_empty(model) -# -# v = MOI.add_variables(model, 2) -# @test MOI.get(model, MOI.NumberOfVariables()) == 2 -# -# cf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0,1.0], v), 0.0) -# c = MOI.add_constraint(model, cf, MOI.LessThan(1.0)) -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}}()) == 1 -# -# vc1 = MOI.add_constraint(model, MOI.SingleVariable(v[1]), MOI.GreaterThan(0.0)) -# # test fallback -# vc2 = MOI.add_constraint(model, v[2], MOI.GreaterThan(0.0)) -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 2 -# -# # note: adding some redundant zero coefficients to catch solvers that don't handle duplicate coefficients correctly: -# objf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([0.0,0.0,-1.0,0.0,0.0,0.0], [v; v; v]), 0.0) -# MOI.set(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), objf) -# MOI.set(model, MOI.ObjectiveSense(), MOI.MinSense) -# -# @test MOI.get(model, MOI.ObjectiveSense()) == MOI.MinSense -# -# if config.query -# vrs = MOI.get(model, MOI.ListOfVariableIndices()) -# @test vrs == v || vrs == reverse(v) -# -# @test objf ≈ MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) -# -# @test cf ≈ MOI.get(model, MOI.ConstraintFunction(), c) -# -# s = MOI.get(model, MOI.ConstraintSet(), c) -# @test s == MOI.LessThan(1.0) -# -# s = MOI.get(model, MOI.ConstraintSet(), vc1) -# @test s == MOI.GreaterThan(0.0) -# -# s = MOI.get(model, MOI.ConstraintSet(), vc2) -# @test s == MOI.GreaterThan(0.0) -# end -# -# if config.solve -# MOI.optimize!(model) -# -# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success -# -# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint -# -# @test MOI.get(model, MOI.ObjectiveValue()) ≈ -1 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [1, 0] atol=atol rtol=rtol -# -# # @test MOI.get(model, MOI.ConstraintPrimal(), c) ≈ 1 atol=atol rtol=rtol -# -# if config.duals -# @test MOI.get(model, MOI.DualStatus()) == MOI.FeasiblePoint -# @test MOI.get(model, MOI.ConstraintDual(), c) ≈ -1 atol=atol rtol=rtol -# -# # reduced costs -# @test MOI.get(model, MOI.ConstraintDual(), vc1) ≈ 0 atol=atol rtol=rtol -# @test MOI.get(model, MOI.ConstraintDual(), vc2) ≈ 1 atol=atol rtol=rtol -# end -# end + @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success + + @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint + + @test MOI.get(model, MOI.ObjectiveValue()) ≈ -1 atol=atol rtol=rtol + + @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [1, 0] atol=atol rtol=rtol + + # @test MOI.get(model, MOI.ConstraintPrimal(), c) ≈ 1 atol=atol rtol=rtol + + if config.duals + @test MOI.get(model, MOI.DualStatus()) == MOI.FeasiblePoint + @test MOI.get(model, MOI.ConstraintDual(), c) ≈ -1 atol=atol rtol=rtol + + # reduced costs + @test MOI.get(model, MOI.ConstraintDual(), vc1) ≈ 0 atol=atol rtol=rtol + @test MOI.get(model, MOI.ConstraintDual(), vc2) ≈ 1 atol=atol rtol=rtol + end +end # # # change objective to Max +x # From 66ee10e2c667e03bbbf550266de197c051dcfeeb Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Sun, 30 Sep 2018 23:42:01 -0400 Subject: [PATCH 09/10] MOI contlinear/contconic tests passing modulo SDP and singularity of A --- Manifest.toml | 2 +- examples/lp/lp.jl | 3 +- examples/namedpoly/namedpoly.jl | 3 +- src/cone.jl | 16 +- src/interpolation.jl | 5 +- src/mathoptinterface.jl | 249 +++++++------ src/nativeinterface.jl | 9 +- src/primitivecones/positivesemidefinite.jl | 12 +- test/runtests.jl | 414 ++------------------- 9 files changed, 200 insertions(+), 513 deletions(-) diff --git a/Manifest.toml b/Manifest.toml index e066ecad1..f5fad534c 100644 --- a/Manifest.toml +++ b/Manifest.toml @@ -141,7 +141,7 @@ deps = ["Distributed", "InteractiveUtils", "Logging", "Random"] uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [[UUIDs]] -deps = ["Random"] +deps = ["Random", "SHA"] uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4" [[Unicode]] diff --git a/examples/lp/lp.jl b/examples/lp/lp.jl index 8593aaaa9..e64fb206d 100644 --- a/examples/lp/lp.jl +++ b/examples/lp/lp.jl @@ -1,5 +1,6 @@ #= -Copyright 2018, David Papp, Sercan Yildiz, and contributors +Copyright 2018, Chris Coey and contributors +Copyright 2018, David Papp, Sercan Yildiz modified from https://github.com/dpapp-github/alfonso/blob/master/random_lp.m solves a simple LP min c'x s.t. Ax = b, x >= 0 diff --git a/examples/namedpoly/namedpoly.jl b/examples/namedpoly/namedpoly.jl index 657e6117d..84dca681b 100644 --- a/examples/namedpoly/namedpoly.jl +++ b/examples/namedpoly/namedpoly.jl @@ -1,5 +1,6 @@ #= -Copyright 2018, David Papp, Sercan Yildiz, and contributors +Copyright 2018, Chris Coey and contributors +Copyright 2018, David Papp, Sercan Yildiz modified from https://github.com/dpapp-github/alfonso/blob/master/polyOpt.m formulates and solves the polynomial optimization problem for a given polynomial, described in the paper: diff --git a/src/cone.jl b/src/cone.jl index 295c09631..d98d435f6 100644 --- a/src/cone.jl +++ b/src/cone.jl @@ -1,25 +1,11 @@ -#= -cone object -=# - +# cone object abstract type PrimitiveCone end # TODO reorder primitive cones so easiest ones to check incone are first mutable struct Cone prms::Vector{PrimitiveCone} idxs::Vector{AbstractVector{Int}} - - function Cone(prms::Vector{PrimitiveCone}, idxs::Vector{AbstractVector{Int}}) - @assert length(prms) == length(idxs) - for k in eachindex(prms) - @assert dimension(prms[k]) == length(idxs[k]) - end - cone = new() - cone.prms = prms - cone.idxs = idxs - return cone - end end Cone() = Cone(PrimitiveCone[], AbstractVector{Int}[]) diff --git a/src/interpolation.jl b/src/interpolation.jl index a9c4cd8ff..bab8c96b3 100644 --- a/src/interpolation.jl +++ b/src/interpolation.jl @@ -1,5 +1,6 @@ #= -Copyright 2018, David Papp, Sercan Yildiz, and contributors +Copyright 2018, Chris Coey and contributors +Copyright 2018, David Papp, Sercan Yildiz modified from files in https://github.com/dpapp-github/alfonso/ https://github.com/dpapp-github/alfonso/blob/master/ChebInterval.m @@ -182,7 +183,7 @@ function approxfekete_data(n::Int, d::Int, calc_w::Bool) pts = ipts[keep_pnt,:] # subset of points indexed with the support of w P0 = M[keep_pnt,1:L] # subset of polynomial evaluations up to total degree d - P = Array(qr(P0).Q) + P = Array(qr(P0).Q) if calc_w w = UpperTriangular(F.R[:,1:U])\(F.Q'*m) diff --git a/src/mathoptinterface.jl b/src/mathoptinterface.jl index 342f79d9b..5758e8c55 100644 --- a/src/mathoptinterface.jl +++ b/src/mathoptinterface.jl @@ -1,19 +1,27 @@ mutable struct Optimizer <: MOI.AbstractOptimizer alf::AlfonsoOpt - + c::Vector{Float64} # linear cost vector, size n + A::AbstractMatrix{Float64} # equality constraint matrix, size p*n + b::Vector{Float64} # equality constraint vector, size p + G::AbstractMatrix{Float64} # cone constraint matrix, size q*n + h::Vector{Float64} # cone constraint vector, size q + cone::Cone # primal constraint cone object objsense::MOI.OptimizationSense - constroffsets::Vector{Int} + objconst::Float64 numeqconstrs::Int - - x - s - y - z - status - solvetime - pobj - dobj + constroffseteq::Vector{Int} + constrprimeq::Vector{Float64} + constroffsetcone::Vector{Int} + constrprimcone::Vector{Float64} + x::Vector{Float64} + s::Vector{Float64} + y::Vector{Float64} + z::Vector{Float64} + status::Symbol + solvetime::Float64 + pobj::Float64 + dobj::Float64 function Optimizer(alf::AlfonsoOpt) opt = new() @@ -62,16 +70,12 @@ conefrommoi(s::MOI.PositiveSemidefiniteConeTriangle) = PositiveSemidefiniteCone( conefrommoi(s::MOI.ExponentialCone) = ExponentialCone() # conefrommoi(s::MOI.PowerCone) = PowerCone(s.exponent) - # build representation as min c'x s.t. Ax = b, x in K # TODO what if some variables x are in multiple cone constraints? function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_attributes=true) @assert !copy_names idxmap = Dict{MOI.Index, MOI.Index}() - # model - # mattr_src = MOI.get(src, MOI.ListOfModelAttributesSet()) - # variables n = MOI.get(src, MOI.NumberOfVariables()) # columns of A # vattr_src = MOI.get(src, MOI.ListOfVariableAttributesSet()) @@ -83,63 +87,68 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ @assert j == n # objective function + obj = MOI.get(src, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) (Jc, Vc) = (Int[], Float64[]) - for t in (MOI.get(src, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}())).terms + for t in obj.terms push!(Jc, idxmap[t.variable_index].value) push!(Vc, t.coefficient) end - ismax = (MOI.get(src, MOI.ObjectiveSense()) == MOI.MaxSense) ? true : false - if ismax - Vc .*= -1 + if MOI.get(src, MOI.ObjectiveSense()) == MOI.MaxSense + Vc .*= -1.0 end + opt.objconst = obj.constant + opt.objsense = MOI.get(src, MOI.ObjectiveSense()) + opt.c = Vector(sparsevec(Jc, Vc, n)) + # constraints getsrccons(F, S) = MOI.get(src, MOI.ListOfConstraintIndices{F,S}()) getconfun(conidx) = MOI.get(src, MOI.ConstraintFunction(), conidx) getconset(conidx) = MOI.get(src, MOI.ConstraintSet(), conidx) - - # pass over vector constraints to construct conic model data - constroffsets = Vector{Int}() - i = 0 # MOI constraint objects + i = 0 # MOI constraint index # equality constraints p = 0 # rows of A (equality constraint matrix) (IA, JA, VA) = (Int[], Int[], Float64[]) (Ib, Vb) = (Int[], Float64[]) + (Icpe, Vcpe) = (Int[], Float64[]) # constraint set constants for opt.constrprimeq + constroffseteq = Vector{Int}() for S in (MOI.EqualTo{Float64}, MOI.Zeros) # TODO can preprocess variables equal to constant for ci in getsrccons(MOI.SingleVariable, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, S}(i) - push!(constroffsets, p) - + push!(constroffseteq, p) p += 1 push!(IA, p) push!(JA, idxmap[getconfun(ci).variable].value) - push!(VA, 1.0) + push!(VA, -1.0) if S == MOI.EqualTo{Float64} push!(Ib, p) - push!(Vb, getconset(ci).value) + push!(Vb, -getconset(ci).value) + push!(Icpe, p) + push!(Vcpe, getconset(ci).value) end end for ci in getsrccons(MOI.ScalarAffineFunction{Float64}, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, S}(i) - push!(constroffsets, p) - + push!(constroffseteq, p) p += 1 fi = getconfun(ci) for vt in fi.terms push!(IA, p) push!(JA, idxmap[vt.variable_index].value) - push!(VA, vt.coefficient) + push!(VA, -vt.coefficient) end push!(Ib, p) if S == MOI.EqualTo{Float64} - push!(Vb, getconset(ci).value - fi.constant) + push!(Vb, fi.constant - getconset(ci).value) + push!(Icpe, p) + push!(Vcpe, getconset(ci).value) else - push!(Vb, -fi.constant) + push!(Vb, fi.constant) end end end @@ -148,38 +157,44 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.VectorOfVariables, MOI.Zeros) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, MOI.Zeros}(i) - push!(constroffsets, p) - + push!(constroffseteq, p) fi = getconfun(ci) dim = MOI.output_dimension(fi) append!(IA, p+1:p+dim) append!(JA, idxmap[vi].value for vi in fi.variables) - append!(VA, ones(dim)) + append!(VA, -ones(dim)) end for ci in getsrccons(MOI.VectorAffineFunction{Float64}, MOI.Zeros) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, MOI.Zeros}(i) - push!(constroffsets, p) - + push!(constroffseteq, p) fi = getconfun(ci) dim = MOI.output_dimension(fi) for vt in fi.terms push!(IA, p + vt.output_index) push!(JA, idxmap[vt.scalar_term.variable_index].value) - push!(VA, vt.scalar_term.coefficient) + push!(VA, -vt.scalar_term.coefficient) end append!(Ib, p+1:p+dim) - append!(Vb, -fi.constants) + append!(Vb, fi.constants) p += dim end - numeqconstrs = i + push!(constroffseteq, p) + opt.A = Matrix(sparse(IA, JA, VA, p, n)) + # opt.A = dropzeros!(sparse(IA, JA, VA, p, n)) + opt.b = Vector(sparsevec(Ib, Vb, p)) + opt.numeqconstrs = i + opt.constrprimeq = Vector(sparsevec(Icpe, Vcpe, p)) + opt.constroffseteq = constroffseteq # conic constraints q = 0 # rows of G (cone constraint matrix) (IG, JG, VG) = (Int[], Int[], Float64[]) (Ih, Vh) = (Int[], Float64[]) + (Icpc, Vcpc) = (Int[], Float64[]) # constraint set constants for opt.constrprimeq + constroffsetcone = Vector{Int}() cone = Cone() # LP constraints: build up one nonnegative cone and/or one nonpositive cone @@ -192,8 +207,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.SingleVariable, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.SingleVariable, S}(i) - push!(constroffsets, q) - + push!(constroffsetcone, q) q += 1 push!(IG, q) push!(JG, idxmap[getconfun(ci).variable].value) @@ -202,17 +216,19 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ if S == MOI.GreaterThan{Float64} push!(Vh, -getconset(ci).lower) push!(nonnegrows, q) + push!(Vcpc, getconset(ci).lower) else push!(Vh, -getconset(ci).upper) push!(nonposrows, q) + push!(Vcpc, getconset(ci).upper) end + push!(Icpc, q) end for ci in getsrccons(MOI.ScalarAffineFunction{Float64}, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, S}(i) - push!(constroffsets, q) - + push!(constroffsetcone, q) q += 1 fi = getconfun(ci) for vt in fi.terms @@ -222,12 +238,15 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ end push!(Ih, q) if S == MOI.GreaterThan{Float64} - push!(Vh, -getconset(ci).lower + fi.constant) + push!(Vh, fi.constant - getconset(ci).lower) push!(nonnegrows, q) + push!(Vcpc, getconset(ci).lower) else - push!(Vh, -getconset(ci).upper + fi.constant) + push!(Vh, fi.constant - getconset(ci).upper) push!(nonposrows, q) + push!(Vcpc, getconset(ci).upper) end + push!(Icpc, q) end end @@ -235,8 +254,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.VectorOfVariables, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, S}(i) - push!(constroffsets, q) - + push!(constroffsetcone, q) for vj in getconfun(ci).variables q += 1 push!(IG, q) @@ -253,8 +271,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.VectorAffineFunction{Float64}, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, S}(i) - push!(constroffsets, q) - + push!(constroffsetcone, q) fi = getconfun(ci) dim = MOI.output_dimension(fi) for vt in fi.terms @@ -282,8 +299,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.VectorOfVariables, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorOfVariables, S}(i) - push!(constroffsets, q) - + push!(constroffsetcone, q) fi = getconfun(ci) dim = MOI.output_dimension(fi) append!(IG, q+1:q+dim) @@ -296,8 +312,7 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ for ci in getsrccons(MOI.VectorAffineFunction{Float64}, S) i += 1 idxmap[ci] = MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64}, S}(i) - push!(constroffsets, q) - + push!(constroffsetcone, q) fi = getconfun(ci) dim = MOI.output_dimension(fi) for vt in fi.terms @@ -312,35 +327,32 @@ function MOI.copy_to(opt::Optimizer, src::MOI.ModelLike; copy_names=false, warn_ end end - # finalize data and load into Alfonso model - c = Vector(dropzeros!(sparsevec(Jc, Vc, n))) - A = dropzeros!(sparse(IA, JA, VA, p, n)) - b = Vector(dropzeros!(sparsevec(Ib, Vb, p))) - G = dropzeros!(sparse(IG, JG, VG, q, n)) - h = Vector(dropzeros!(sparsevec(Ih, Vh, q))) - # TODO make converting to dense optional - # Alfonso.load_data!(opt.alf, c, A, b, G, h, cone) # sparse - Alfonso.load_data!(opt.alf, c, Matrix(A), b, Matrix(G), h, cone) # dense - - # store information needed by MOI getter functions - opt.objsense = MOI.get(src, MOI.ObjectiveSense()) - opt.constroffsets = constroffsets - opt.numeqconstrs = numeqconstrs + push!(constroffsetcone, q) + opt.G = Matrix(sparse(IG, JG, VG, q, n)) + # opt.G = dropzeros!(sparse(IG, JG, VG, q, n)) + opt.h = Vector(sparsevec(Ih, Vh, q)) + opt.cone = cone + opt.constroffsetcone = constroffsetcone + opt.constrprimcone = Vector(sparsevec(Icpc, Vcpc, q)) return idxmap end function MOI.optimize!(opt::Optimizer) - solve!(opt.alf) + Alfonso.load_data!(opt.alf, opt.c, opt.A, opt.b, opt.G, opt.h, opt.cone) # dense + Alfonso.solve!(opt.alf) - opt.x = get_x(opt.alf) - opt.s = get_s(opt.alf) - opt.y = get_y(opt.alf) - opt.z = get_z(opt.alf) - opt.status = get_status(opt.alf) - opt.solvetime = get_solvetime(opt.alf) - opt.pobj = get_pobj(opt.alf) - opt.dobj = get_dobj(opt.alf) + opt.x = Alfonso.get_x(opt.alf) + opt.constrprimeq += opt.b - opt.A*opt.x + opt.s = Alfonso.get_s(opt.alf) + opt.constrprimcone += opt.s + opt.y = Alfonso.get_y(opt.alf) + opt.z = Alfonso.get_z(opt.alf) + + opt.status = Alfonso.get_status(opt.alf) + opt.solvetime = Alfonso.get_solvetime(opt.alf) + opt.pobj = Alfonso.get_pobj(opt.alf) + opt.dobj = Alfonso.get_dobj(opt.alf) return nothing end @@ -362,9 +374,9 @@ end function MOI.get(opt::Optimizer, ::MOI.ObjectiveValue) if opt.objsense == MOI.MinSense - return opt.pobj + return opt.pobj + opt.objconst elseif opt.objsense == MOI.MaxSense - return -opt.pobj + return -opt.pobj + opt.objconst else error("no objective sense is set") end @@ -372,9 +384,9 @@ end function MOI.get(opt::Optimizer, ::MOI.ObjectiveBound) if opt.objsense == MOI.MinSense - return opt.dobj + return opt.dobj + opt.objconst elseif opt.objsense == MOI.MaxSense - return -opt.dobj + return -opt.dobj + opt.objconst else error("no objective sense is set") end @@ -418,47 +430,58 @@ end MOI.get(opt::Optimizer, ::MOI.VariablePrimal, vi::MOI.VariableIndex) = opt.x[vi.value] MOI.get(opt::Optimizer, a::MOI.VariablePrimal, vi::Vector{MOI.VariableIndex}) = MOI.get.(opt, a, vi) -function MOI.get(opt::Optimizer, ::MOI.ConstraintDual, ci::MOI.ConstraintIndex{F, S}) where {F, S} - os = opt.constroffsets +function MOI.get(opt::Optimizer, ::MOI.ConstraintDual, ci::MOI.ConstraintIndex{F, S}) where {F <: MOI.AbstractFunction, S <: MOI.AbstractScalarSet} + # scalar set + i = ci.value + if i <= opt.numeqconstrs + # constraint is an equality + return opt.y[opt.constroffseteq[i]+1] + else + # constraint is conic + i -= opt.numeqconstrs + return opt.z[opt.constroffsetcone[i]+1] + end +end +function MOI.get(opt::Optimizer, ::MOI.ConstraintDual, ci::MOI.ConstraintIndex{F, S}) where {F <: MOI.AbstractFunction, S <: MOI.AbstractVectorSet} + # vector set i = ci.value if i <= opt.numeqconstrs # constraint is an equality + os = opt.constroffseteq return opt.y[os[i]+1:os[i+1]] else # constraint is conic + i -= opt.numeqconstrs + os = opt.constroffsetcone return opt.z[os[i]+1:os[i+1]] end end MOI.get(opt::Optimizer, a::MOI.ConstraintDual, ci::Vector{MOI.ConstraintIndex}) = MOI.get.(opt, a, ci) -# TODO use x and s. if conic then just return s. otherwise have to use x - rhs after saving rhs in opt maybe -MOI.get(opt::Optimizer, ::MOI.ConstraintPrimal, ci::MOI.ConstraintIndex) = NaN +function MOI.get(opt::Optimizer, ::MOI.ConstraintPrimal, ci::MOI.ConstraintIndex{F, S}) where {F <: MOI.AbstractFunction, S <: MOI.AbstractScalarSet} + # scalar set + i = ci.value + if i <= opt.numeqconstrs + # constraint is an equality + return opt.constrprimeq[opt.constroffseteq[i]+1] + else + # constraint is conic + i -= opt.numeqconstrs + return opt.constrprimcone[opt.constroffsetcone[i]+1] + end +end +function MOI.get(opt::Optimizer, ::MOI.ConstraintPrimal, ci::MOI.ConstraintIndex{F, S}) where {F <: MOI.AbstractFunction, S <: MOI.AbstractVectorSet} + # vector set + i = ci.value + if i <= opt.numeqconstrs + # constraint is an equality + os = opt.constroffseteq + return opt.constrprimeq[os[i]+1:os[i+1]] + else + # constraint is conic + i -= opt.numeqconstrs + os = opt.constroffsetcone + return opt.constrprimcone[os[i]+1:os[i+1]] + end +end MOI.get(opt::Optimizer, a::MOI.ConstraintPrimal, ci::Vector{MOI.ConstraintIndex}) = MOI.get.(opt, a, ci) - - - - - -# function MOI.get(instance::Optimizer, ::MOI.ConstraintPrimal, ci::CI{<:MOI.AbstractFunction, MOI.Zeros}) -# rows = constrrows(instance, ci) -# zeros(length(rows)) -# end -# function MOI.get(instance::Optimizer, ::MOI.ConstraintPrimal, ci::CI{<:MOI.AbstractFunction, <:MOI.EqualTo}) -# offset = constroffset(instance, ci) -# setconstant(instance, offset, ci) -# end -# function MOI.get(instance::Optimizer, ::MOI.ConstraintPrimal, ci::CI{<:MOI.AbstractFunction, S}) where S <: MOI.AbstractSet -# offset = constroffset(instance, ci) -# rows = constrrows(instance, ci) -# _unshift(instance, offset, scalecoef(rows, reorderval(instance.sol.slack[offset .+ rows], S), false, S), ci) -# end -# -# _dual(instance, ci::CI{<:MOI.AbstractFunction, <:ZeroCones}) = instance.sol.dual_eq -# _dual(instance, ci::CI) = instance.sol.dual_ineq -# function MOI.get(instance::Optimizer, ::MOI.ConstraintDual, ci::CI{<:MOI.AbstractFunction, S}) where S <: MOI.AbstractSet -# offset = constroffset(instance, ci) -# rows = constrrows(instance, ci) -# scalecoef(rows, reorderval(_dual(instance, ci)[offset .+ rows], S), false, S) -# end -# -# MOI.get(instance::Optimizer, ::MOI.ResultCount) = 1 diff --git a/src/nativeinterface.jl b/src/nativeinterface.jl index d3c24c55a..b14c98318 100644 --- a/src/nativeinterface.jl +++ b/src/nativeinterface.jl @@ -191,6 +191,11 @@ function load_data!( error("number of constraint rows is not consistent in G and h") end + @assert length(cone.prms) == length(cone.idxs) + for k in eachindex(cone.prms) + @assert dimension(cone.prms[k]) == length(cone.idxs[k]) + end + # perform QR decomposition of A' for use in linear system solves # TODO reduce allocs, improve efficiency # A' = [Q1 Q2] * [R1; 0] @@ -574,7 +579,7 @@ end function findinitialiterate!(tx::Vector{Float64}, ty::Vector{Float64}, tz::Vector{Float64}, ts::Vector{Float64}, ls_ts::Vector{Float64}, bnu::Float64, alf::AlfonsoOpt) (b, G, h, cone) = (alf.b, alf.G, alf.h, alf.cone) L = alf.L - (Q2, RiQ1, GQ2, Q2GHGQ2, bxGHbz, Q1x, rhs, Q2div, Q2x, HGxi) = (L.Q2, L.RiQ1, L.GQ2, L.Q2GHGQ2, L.bxGHbz, L.Q1x, L.rhs, L.Q2div, L.Q2x, L.HGxi) + (Q2, RiQ1, GQ2, GHGQ2, Q2GHGQ2, bxGHbz, Q1x, rhs, Q2div, Q2x, HGxi) = (L.Q2, L.RiQ1, L.GQ2, L.GHGQ2, L.Q2GHGQ2, L.bxGHbz, L.Q1x, L.rhs, L.Q2div, L.Q2x, L.HGxi) alf.verbose && println("\nfinding initial iterate") @@ -642,7 +647,7 @@ function findinitialiterate!(tx::Vector{Float64}, ty::Vector{Float64}, tz::Vecto mu = (dot(tz, ts) + tau*kap)/bnu @assert abs(1.0 - mu) < 1e-10 # @assert calcnbhd(tau*kap, mu, copy(tz), copy(tz), cone) < 1e-6 - + alf.verbose && println("initial iterate found") return (tau, kap, mu) diff --git a/src/primitivecones/positivesemidefinite.jl b/src/primitivecones/positivesemidefinite.jl index 212d5c635..b66902b7d 100644 --- a/src/primitivecones/positivesemidefinite.jl +++ b/src/primitivecones/positivesemidefinite.jl @@ -25,7 +25,17 @@ end dimension(prm::PositiveSemidefiniteCone) = prm.dim barrierpar_prm(prm::PositiveSemidefiniteCone) = prm.side -getintdir_prm!(arr::AbstractVector{Float64}, prm::PositiveSemidefiniteCone) = mattovec!(arr, Matrix(1.0I, prm.side, prm.side)) # TODO eliminate allocs +function getintdir_prm!(arr::AbstractVector{Float64}, prm::PositiveSemidefiniteCone) + for i in 1:prm.side, j in i:prm.side + if i == j + prm.mat[i,j] = 1.0 + else + prm.mat[i,j] = prm.mat[j,i] = 0.0 + end + end + mattovec!(arr, prm.mat) + return arr +end loadpnt_prm!(prm::PositiveSemidefiniteCone, pnt::AbstractVector{Float64}) = (prm.pnt = pnt) function incone_prm(prm::PositiveSemidefiniteCone) diff --git a/test/runtests.jl b/test/runtests.jl index 7d42971c5..22240e540 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2,19 +2,24 @@ using Alfonso using Test -# verbflag = true # Alfonso verbose option -# -# # TODO interpolation tests -# -# # load optimizer builder functions from examples folder -# egs_dir = joinpath(@__DIR__, "../examples") -# include(joinpath(egs_dir, "envelope/envelope.jl")) -# include(joinpath(egs_dir, "lp/lp.jl")) -# include(joinpath(egs_dir, "namedpoly/namedpoly.jl")) -# -# # run native and MOI interfaces on examples -# include(joinpath(@__DIR__, "nativeexamples.jl")) +# native interface tests + +verbflag = false # Alfonso verbose option + +# TODO interpolation tests + +# load optimizer builder functions from examples folder +egs_dir = joinpath(@__DIR__, "../examples") +include(joinpath(egs_dir, "envelope/envelope.jl")) +include(joinpath(egs_dir, "lp/lp.jl")) +include(joinpath(egs_dir, "namedpoly/namedpoly.jl")) + +# run native and MOI interfaces on examples +include(joinpath(@__DIR__, "nativeexamples.jl")) + + +# MathOptInterface tests import MathOptInterface MOI = MathOptInterface @@ -44,9 +49,10 @@ MOIU.@model(AlfonsoModelData, optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Optimizer()) +# TODO tols relaxed temporarily from 1e-4 config = MOIT.TestConfig( - atol=1e-4, - rtol=1e-4, + atol=1e-3, + rtol=1e-3, solve=true, query=true, modify_lhs=true, @@ -54,371 +60,25 @@ config = MOIT.TestConfig( infeas_certificates=true, ) -# @testset "Continuous linear problems" begin -# MOIT.contlineartest( -# MOIB.SplitInterval{Float64}( -# optimizer -# ), -# config) -# end - -# @testset "Continuous conic problems" begin -# exclude = ["rootdet", "logdet"] # TODO bridges not working? should not need to exclude in future -# MOIT.contconictest( -# MOIB.SquarePSD{Float64}( -# MOIB.GeoMean{Float64}( -# # MOIB.LogDet{Float64}( -# # MOIB.RootDet{Float64}( -# optimizer -# )),#)), -# config, exclude) -# end - - - -model = optimizer - -atol = config.atol -rtol = config.rtol - -# simple 2 variable, 1 constraint problem -# min -x -# st x + y <= 1 (x + y - 1 ∈ Nonpositives) -# x, y >= 0 (x, y ∈ Nonnegatives) - -@test MOI.supports(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) -@test MOI.supports(model, MOI.ObjectiveSense()) -@test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}) -@test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.GreaterThan{Float64}) -@test MOI.supports_constraint(model, MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}) -@test MOI.supports_constraint(model, MOI.SingleVariable, MOI.EqualTo{Float64}) -@test MOI.supports_constraint(model, MOI.SingleVariable, MOI.GreaterThan{Float64}) - -MOI.empty!(model) -@test MOI.is_empty(model) - -v = MOI.add_variables(model, 2) -@test MOI.get(model, MOI.NumberOfVariables()) == 2 - -cf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0,1.0], v), 0.0) -c = MOI.add_constraint(model, cf, MOI.LessThan(1.0)) -@test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}}()) == 1 - -vc1 = MOI.add_constraint(model, MOI.SingleVariable(v[1]), MOI.GreaterThan(0.0)) -# test fallback -vc2 = MOI.add_constraint(model, v[2], MOI.GreaterThan(0.0)) -@test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 2 - -# note: adding some redundant zero coefficients to catch solvers that don't handle duplicate coefficients correctly: -objf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([0.0,0.0,-1.0,0.0,0.0,0.0], [v; v; v]), 0.0) -MOI.set(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), objf) -MOI.set(model, MOI.ObjectiveSense(), MOI.MinSense) - -@test MOI.get(model, MOI.ObjectiveSense()) == MOI.MinSense - -if config.query - vrs = MOI.get(model, MOI.ListOfVariableIndices()) - @test vrs == v || vrs == reverse(v) - - @test objf ≈ MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) - - @test cf ≈ MOI.get(model, MOI.ConstraintFunction(), c) - - s = MOI.get(model, MOI.ConstraintSet(), c) - @test s == MOI.LessThan(1.0) - - s = MOI.get(model, MOI.ConstraintSet(), vc1) - @test s == MOI.GreaterThan(0.0) - - s = MOI.get(model, MOI.ConstraintSet(), vc2) - @test s == MOI.GreaterThan(0.0) +@testset "Continuous linear problems" begin + MOIT.contlineartest( + MOIB.SplitInterval{Float64}( + optimizer + ), + config) end -if config.solve - MOI.optimize!(model) - - @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success - - @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint - - @test MOI.get(model, MOI.ObjectiveValue()) ≈ -1 atol=atol rtol=rtol - - @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [1, 0] atol=atol rtol=rtol - - # @test MOI.get(model, MOI.ConstraintPrimal(), c) ≈ 1 atol=atol rtol=rtol - - if config.duals - @test MOI.get(model, MOI.DualStatus()) == MOI.FeasiblePoint - @test MOI.get(model, MOI.ConstraintDual(), c) ≈ -1 atol=atol rtol=rtol - - # reduced costs - @test MOI.get(model, MOI.ConstraintDual(), vc1) ≈ 0 atol=atol rtol=rtol - @test MOI.get(model, MOI.ConstraintDual(), vc2) ≈ 1 atol=atol rtol=rtol - end +@testset "Continuous conic problems" begin + exclude = ["rootdet", "logdet"] # TODO bridges not working? should not need to exclude in future + MOIT.contconictest( + MOIB.SquarePSD{Float64}( + MOIB.GeoMean{Float64}( + # MOIB.LogDet{Float64}( + # MOIB.RootDet{Float64}( + optimizer + )),#)), + config, exclude) end -# -# # change objective to Max +x -# -# objf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0,0.0], v), 0.0) -# MOI.set(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), objf) -# MOI.set(model, MOI.ObjectiveSense(), MOI.MaxSense) -# -# if config.query -# @test objf ≈ MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) -# end -# -# @test MOI.get(model, MOI.ObjectiveSense()) == MOI.MaxSense -# -# if config.solve -# MOI.optimize!(model) -# -# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success -# -# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint -# -# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 1 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [1, 0] atol=atol rtol=rtol -# -# if config.duals -# @test MOI.get(model, MOI.DualStatus()) == MOI.FeasiblePoint -# @test MOI.get(model, MOI.ConstraintDual(), c) ≈ -1 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.ConstraintDual(), vc1) ≈ 0 atol=atol rtol=rtol -# @test MOI.get(model, MOI.ConstraintDual(), vc2) ≈ 1 atol=atol rtol=rtol -# end -# end -# -# # add new variable to get : -# # max x + 2z -# # s.t. x + y + z <= 1 -# # x,y,z >= 0 -# -# z = MOI.add_variable(model) -# push!(v, z) -# @test v[3] == z -# -# if config.query -# # Test that the modification of v has not affected the model -# vars = map(t -> t.variable_index, MOI.get(model, MOI.ConstraintFunction(), c).terms) -# @test vars == [v[1], v[2]] || vars == [v[2], v[1]] -# @test MOI.ScalarAffineFunction([MOI.ScalarAffineTerm(1.0, v[1])], 0.0) ≈ MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) -# end -# -# vc3 = MOI.add_constraint(model, MOI.SingleVariable(v[3]), MOI.GreaterThan(0.0)) -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 3 -# -# if config.modify_lhs -# MOI.modify(model, c, MOI.ScalarCoefficientChange{Float64}(z, 1.0)) -# else -# MOI.delete(model, c) -# cf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0,1.0,1.0], v), 0.0) -# c = MOI.add_constraint(model, cf, MOI.LessThan(1.0)) -# end -# -# MOI.modify(model, -# MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), -# MOI.ScalarCoefficientChange{Float64}(z, 2.0) -# ) -# -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}}()) == 1 -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 3 -# -# if config.solve -# MOI.optimize!(model) -# -# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success -# -# @test MOI.get(model, MOI.ResultCount()) >= 1 -# -# @test MOI.get(model, MOI.PrimalStatus(1)) == MOI.FeasiblePoint -# -# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 2 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [0, 0, 1] atol=atol rtol=rtol -# -# # @test MOI.get(model, MOI.ConstraintPrimal(), c) ≈ 1 atol=atol rtol=rtol -# -# if config.duals -# @test MOI.get(model, MOI.DualStatus()) == MOI.FeasiblePoint -# @test MOI.get(model, MOI.ConstraintDual(), c) ≈ -2 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.ConstraintDual(), vc1) ≈ 1 atol=atol rtol=rtol -# @test MOI.get(model, MOI.ConstraintDual(), vc2) ≈ 2 atol=atol rtol=rtol -# @test MOI.get(model, MOI.ConstraintDual(), vc3) ≈ 0 atol=atol rtol=rtol -# end -# end -# -# # setting lb of x to -1 to get : -# # max x + 2z -# # s.t. x + y + z <= 1 -# # x >= -1 -# # y,z >= 0 -# MOI.set(model, MOI.ConstraintSet(), vc1, MOI.GreaterThan(-1.0)) -# -# if config.solve -# MOI.optimize!(model) -# -# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success -# -# @test MOI.get(model, MOI.ResultCount()) >= 1 -# -# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint -# -# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 3 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [-1, 0, 2] atol=atol rtol=rtol -# end -# -# # put lb of x back to 0 and fix z to zero to get : -# # max x + 2z -# # s.t. x + y + z <= 1 -# # x, y >= 0, z = 0 (vc3) -# MOI.set(model, MOI.ConstraintSet(), vc1, MOI.GreaterThan(0.0)) -# -# MOI.delete(model, vc3) -# -# vc3 = MOI.add_constraint(model, MOI.SingleVariable(v[3]), MOI.EqualTo(0.0)) -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 2 -# -# if config.solve -# MOI.optimize!(model) -# -# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success -# -# @test MOI.get(model, MOI.ResultCount()) >= 1 -# -# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint -# -# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 1 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [1, 0, 0] atol=atol rtol=rtol -# end -# -# # modify affine linear constraint set to be == 2 to get : -# # max x + 2z -# # s.t. x + y + z == 2 (c) -# # x,y >= 0, z = 0 -# MOI.delete(model, c) -# # note: adding some redundant zero coefficients to catch solvers that don't handle duplicate coefficients correctly: -# cf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([0.0,0.0,0.0,1.0,1.0,1.0,0.0,0.0,0.0], [v; v; v]), 0.0) -# c = MOI.add_constraint(model, cf, MOI.EqualTo(2.0)) -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}}()) == 0 -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.EqualTo{Float64}}()) == 1 -# -# if config.solve -# MOI.optimize!(model) -# -# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success -# -# @test MOI.get(model, MOI.ResultCount()) >= 1 -# -# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint -# -# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 2 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [2, 0, 0] atol=atol rtol=rtol -# end -# -# # modify objective function to x + 2y to get : -# # max x + 2y -# # s.t. x + y + z == 2 (c) -# # x,y >= 0, z = 0 -# -# objf = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0,2.0,0.0], v), 0.0) -# MOI.set(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), objf) -# MOI.set(model, MOI.ObjectiveSense(), MOI.MaxSense) -# -# if config.query -# @test objf ≈ MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) -# end -# -# if config.solve -# MOI.optimize!(model) -# -# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success -# -# @test MOI.get(model, MOI.ResultCount()) >= 1 -# -# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint -# -# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 4 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [0, 2, 0] atol=atol rtol=rtol -# end -# -# # add constraint x - y >= 0 (c2) to get : -# # max x+2y -# # s.t. x + y + z == 2 (c) -# # x - y >= 0 (c2) -# # x,y >= 0 (vc1,vc2), z = 0 (vc3) -# -# cf2 = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([1.0, -1.0, 0.0], v), 0.0) -# c2 = MOI.add_constraint(model, cf2, MOI.GreaterThan(0.0)) -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.EqualTo{Float64}}()) == 1 -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.GreaterThan{Float64}}()) == 1 -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}}()) == 0 -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.EqualTo{Float64}}()) == 1 -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 2 -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.LessThan{Float64}}()) == 0 -# -# if config.solve -# MOI.optimize!(model) -# -# @test MOI.get(model, MOI.TerminationStatus()) == MOI.Success -# -# @test MOI.get(model, MOI.ResultCount()) >= 1 -# -# @test MOI.get(model, MOI.PrimalStatus()) == MOI.FeasiblePoint -# -# @test MOI.get(model, MOI.ObjectiveValue()) ≈ 3 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.VariablePrimal(), v) ≈ [1, 1, 0] atol=atol rtol=rtol -# -# # @test MOI.get(model, MOI.ConstraintPrimal(), c) ≈ 2 atol=atol rtol=rtol -# # -# # @test MOI.get(model, MOI.ConstraintPrimal(), c2) ≈ 0 atol=atol rtol=rtol -# # -# # @test MOI.get(model, MOI.ConstraintPrimal(), vc1) ≈ 1 atol=atol rtol=rtol -# # -# # @test MOI.get(model, MOI.ConstraintPrimal(), vc2) ≈ 1 atol=atol rtol=rtol -# # -# # @test MOI.get(model, MOI.ConstraintPrimal(), vc3) ≈ 0 atol=atol rtol=rtol -# -# if config.duals -# @test MOI.get(model, MOI.DualStatus(1)) == MOI.FeasiblePoint -# -# @test MOI.get(model, MOI.ConstraintDual(), c) ≈ -1.5 atol=atol rtol=rtol -# @test MOI.get(model, MOI.ConstraintDual(), c2) ≈ 0.5 atol=atol rtol=rtol -# -# @test MOI.get(model, MOI.ConstraintDual(), vc1) ≈ 0 atol=atol rtol=rtol -# @test MOI.get(model, MOI.ConstraintDual(), vc2) ≈ 0 atol=atol rtol=rtol -# @test MOI.get(model, MOI.ConstraintDual(), vc3) ≈ 1.5 atol=atol rtol=rtol -# end -# end -# -# if config.query -# @test MOI.get(model, MOI.ConstraintFunction(), c2) ≈ cf2 -# end -# -# # delete variable x to get : -# # max 2y -# # s.t. y + z == 2 -# # - y >= 0 -# # y >= 0, z = 0 -# -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 2 -# MOI.delete(model, v[1]) -# @test MOI.get(model, MOI.NumberOfConstraints{MOI.SingleVariable,MOI.GreaterThan{Float64}}()) == 1 -# -# if config.query -# @test MOI.get(model, MOI.ConstraintFunction(), c2) ≈ MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([-1.0, 0.0], [v[2], z]), 0.0) -# -# vrs = MOI.get(model, MOI.ListOfVariableIndices()) -# @test vrs == [v[2], z] || vrs == [z, v[2]] -# @test MOI.get(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}()) ≈ MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([2.0, 0.0], [v[2], z]), 0.0) -# end -# -# + return nothing From a8bd54ecfb63cfb0a0d8174ae112bbe70366ea34 Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Mon, 1 Oct 2018 15:17:33 -0400 Subject: [PATCH 10/10] exclude SDP tests for now, restore native tests --- src/nativeinterface.jl | 4 ++-- test/moiexamples.jl | 0 test/{nativeexamples.jl => native.jl} | 0 test/runtests.jl | 11 +++++------ 4 files changed, 7 insertions(+), 8 deletions(-) delete mode 100644 test/moiexamples.jl rename test/{nativeexamples.jl => native.jl} (100%) diff --git a/src/nativeinterface.jl b/src/nativeinterface.jl index b14c98318..ad4e9c06d 100644 --- a/src/nativeinterface.jl +++ b/src/nativeinterface.jl @@ -599,7 +599,7 @@ function findinitialiterate!(tx::Vector{Float64}, ty::Vector{Float64}, tz::Vecto mul!(Q2GHGQ2, GQ2', GQ2) F = cholesky!(Symmetric(Q2GHGQ2), check=false) if !issuccess(F) - alf.verbose && println("linear system matrix was nearly not positive definite") + alf.verbose && println("linear system matrix was not positive definite") mul!(Q2GHGQ2, Q2', GHGQ2) F = bunchkaufman!(Symmetric(Q2GHGQ2)) end @@ -679,7 +679,7 @@ function finddirection!(rhs_tx::Vector{Float64}, rhs_ty::Vector{Float64}, rhs_tz # TODO does it matter that F could be either type? F = cholesky!(Symmetric(Q2GHGQ2), check=false) if !issuccess(F) - alf.verbose && println("linear system matrix was nearly not positive definite") + alf.verbose && println("linear system matrix was not positive definite") mul!(Q2GHGQ2, Q2', GHGQ2) F = bunchkaufman!(Symmetric(Q2GHGQ2)) end diff --git a/test/moiexamples.jl b/test/moiexamples.jl deleted file mode 100644 index e69de29bb..000000000 diff --git a/test/nativeexamples.jl b/test/native.jl similarity index 100% rename from test/nativeexamples.jl rename to test/native.jl diff --git a/test/runtests.jl b/test/runtests.jl index 22240e540..734519c78 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -49,10 +49,9 @@ MOIU.@model(AlfonsoModelData, optimizer = MOIU.CachingOptimizer(AlfonsoModelData{Float64}(), Alfonso.Optimizer()) -# TODO tols relaxed temporarily from 1e-4 config = MOIT.TestConfig( - atol=1e-3, - rtol=1e-3, + atol=1e-4, + rtol=1e-4, solve=true, query=true, modify_lhs=true, @@ -69,14 +68,14 @@ config = MOIT.TestConfig( end @testset "Continuous conic problems" begin - exclude = ["rootdet", "logdet"] # TODO bridges not working? should not need to exclude in future + exclude = ["rootdet", "logdet", "sdp"] # TODO bridges not working? should not need to exclude in future MOIT.contconictest( - MOIB.SquarePSD{Float64}( + # MOIB.SquarePSD{Float64}( MOIB.GeoMean{Float64}( # MOIB.LogDet{Float64}( # MOIB.RootDet{Float64}( optimizer - )),#)), + ),#))), config, exclude) end