Skip to content

Commit

Permalink
Cleanup.
Browse files Browse the repository at this point in the history
  • Loading branch information
tkoolen committed Feb 14, 2018
1 parent 8811208 commit 44af411
Showing 1 changed file with 129 additions and 102 deletions.
231 changes: 129 additions & 102 deletions src/MathOptInterfaceOSQP.jl
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,8 @@ const MOIU = MathOptInterfaceUtilities
const CI = MOI.ConstraintIndex
const VI = MOI.VariableIndex

const SparseTriplets = Tuple{Vector{<:Integer}, Vector{<:Integer}, Vector{<:Any}}

const SingleVariable = MOI.SingleVariable
const Affine = MOI.ScalarAffineFunction{Float64}
const Quadratic = MOI.ScalarQuadraticFunction{Float64}
Expand All @@ -23,7 +25,6 @@ const EqualTo = MOI.EqualTo{Float64}
const IntervalConvertible = Union{Interval, LessThan, GreaterThan, EqualTo}

import OSQP
import MathOptInterfaceUtilities: IndexMap

# TODO: consider moving to MOI:
constant(f::SingleVariable) = 0
Expand Down Expand Up @@ -57,69 +58,18 @@ MOI.isempty(optimizer::OSQPOptimizer) = optimizer.isempty
function MOI.copy!(dest::OSQPOptimizer, src::MOI.ModelLike)
# TODO: cangets (awaiting https://github.com/JuliaOpt/MathOptInterface.jl/pull/210)

# Empty
MOI.empty!(dest)

# Compute problem size, set up index map, and check constraint types
n = MOI.get(src, MOI.NumberOfVariables())
idxmap = IndexMap()
vis_src = MOI.get(src, MOI.ListOfVariableIndices())
for i in eachindex(vis_src)
idxmap[vis_src[i]] = VI(i)
end
m = 0
let i = 0
for (F, S) in MOI.get(src, MOI.ListOfConstraints())
MOI.supportsconstraint(dest, F, S) || return MOI.CopyResult(MOI.CopyUnsupportedConstraint, "Unsupported $F-in-$S constraint", idxmap)
m += MOI.get(src, MOI.NumberOfConstraints{F, S}())
cis_src = MOI.get(src, MOI.ListOfConstraintIndices{F, S}())
for ci in cis_src
i += 1
idxmap[ci] = CI{F, S}(i)
end
end
end

# Allocate storage for problem data
P = spzeros(n, n)
q = zeros(n)
A = spzeros(m, n)
l = zeros(m)
u = zeros(m)

# Process objective function
# prefer the simplest form
sense = dest.sense = MOI.get(src, MOI.ObjectiveSense())
if sense != MOI.FeasibilitySense
if MOI.canget(src, MOI.ObjectiveFunction{SingleVariable}())
objfun = MOI.get(src, MOI.ObjectiveFunction{SingleVariable}())
processobjective!(P, q, objfun, idxmap)
elseif MOI.canget(src, MOI.ObjectiveFunction{Affine}())
objfun = MOI.get(src, MOI.ObjectiveFunction{Affine}())
processobjective!(P, q, objfun, idxmap)
dest.objconstant = constant(objfun)
elseif MOI.canget(src, MOI.ObjectiveFunction{Quadratic}())
objfun = MOI.get(src, MOI.ObjectiveFunction{Quadratic}())
processobjective!(P, q, objfun, idxmap)
dest.objconstant = constant(objfun)
else
return MOI.CopyResult(MOI.CopyOtherError, "No suitable objective function found", idxmap)
end
sense == MOI.MaxSense && (scale!(P, -1); scale!(q, -1); dest.objconstant = -dest.objconstant)
end

# Process constraints
idxmap = MOIU.IndexMap(dest, src)
for (F, S) in MOI.get(src, MOI.ListOfConstraints())
processconstraints!(A, l, u, src, idxmap, F, S)
MOI.supportsconstraint(dest, F, S) || return MOI.CopyResult(MOI.CopyUnsupportedConstraint, "Unsupported $F-in-$S constraint", idxmap)
end

# Load data into OSQP Model
dest.sense, P, q, dest.objconstant = processobjective(src, idxmap)
A, l, u = processconstraints(src, idxmap)
OSQP.setup!(dest.inner; P = P, q = q, A = A, l = l, u = u, dest.settings...)

# Process optimizer attributes
# TODO

# TODO: clean up:
# Process variable attributes
m, n = size(A)
if MOI.canget(src, MOI.VariablePrimalStart(), VI)
x = zeros(n)
i = 1
Expand All @@ -145,67 +95,143 @@ function MOI.copy!(dest::OSQPOptimizer, src::MOI.ModelLike)
return MOI.CopyResult(MOI.CopySuccess, "", idxmap)
end

# These could probably move to MOIU:
function process_affine_objective!(P::SparseMatrixCSC, q::Vector, variables::Vector{VI}, coefficients::Vector, idxmap)
q[:] = 0
n = length(coefficients)
@assert length(variables) == n
for i = 1 : n
q[idxmap[variables[i]].value] += coefficients[i]
"""
Set up index map from `src` variables and constraints to `dest` variables and constraints.
"""
function MOIU.IndexMap(dest::OSQPOptimizer, src::MOI.ModelLike)
idxmap = MOIU.IndexMap()
vis_src = MOI.get(src, MOI.ListOfVariableIndices())
for i in eachindex(vis_src)
idxmap[vis_src[i]] = VI(i)
end
i = 0
for (F, S) in MOI.get(src, MOI.ListOfConstraints())
cis_src = MOI.get(src, MOI.ListOfConstraintIndices{F, S}())
for ci in cis_src
i += 1
idxmap[ci] = CI{F, S}(i)
end
end
idxmap
end

function processobjective!(P::SparseMatrixCSC, q::Vector, objfun::MOI.ScalarAffineFunction, idxmap)
process_affine_objective!(P, q, objfun.variables, objfun.coefficients, idxmap)
end
"""
Return objective sense, as well as matrix `P`, vector `q`, and scalar `c` such that objective function is 1/2 `x' P x + q' x + c`.
"""
function processobjective(src::MOI.ModelLike, idxmap)
sense = MOI.get(src, MOI.ObjectiveSense())
n = MOI.get(src, MOI.NumberOfVariables())

function processobjective!(P::SparseMatrixCSC, q::Vector, objfun::MOI.ScalarQuadraticFunction, idxmap)
process_affine_objective!(P, q, objfun.affine_variables, objfun.affine_coefficients, idxmap)
nquadratic = length(objfun.quadratic_coefficients)
@assert length(objfun.quadratic_rowvariables) == length(objfun.quadratic_colvariables) == nquadratic
P[:] = 0
for i = 1 : nquadratic
row = idxmap[objfun.quadratic_rowvariables[i]].value
col = idxmap[objfun.quadratic_colvariables[i]].value
P[row, col] += objfun.quadratic_coefficients[i]
P[col, row] = P[row, col]
if sense != MOI.FeasibilitySense
if MOI.canget(src, MOI.ObjectiveFunction{SingleVariable}())
fsingle = MOI.get(src, MOI.ObjectiveFunction{SingleVariable}())
P = spzeros(n, n)
q = zeros(n)
q[idxmap[fsingle.variable.value]] = 1
c = 0.
elseif MOI.canget(src, MOI.ObjectiveFunction{Affine}())
faffine = MOI.get(src, MOI.ObjectiveFunction{Affine}())
P = spzeros(n, n)
q = processlinearterm(faffine.variables, faffine.coefficients, idxmap)
c = faffine.constant
elseif MOI.canget(src, MOI.ObjectiveFunction{Quadratic}())
fquadratic = MOI.get(src, MOI.ObjectiveFunction{Quadratic}())
I = [idxmap[var].value for var in fquadratic.quadratic_rowvariables]
J = [idxmap[var].value for var in fquadratic.quadratic_colvariables]
V = fquadratic.quadratic_coefficients
symmetrize!((I, J, V))
P = sparse(I, J, V, n, n)
q = processlinearterm(fquadratic.affine_variables, fquadratic.affine_coefficients, idxmap)
c = fquadratic.constant
else
error("No suitable objective function found")
end
sense == MOI.MaxSense && (scale!(P, -1); scale!(q, -1); c = -c)
else
P = spzeros(n, n)
q = zeros(n)
c = 0.
end
sense, P, q, c
end

function processobjective!(P::SparseMatrixCSC, q::Vector, objfun::MOI.SingleVariable, idxmap)
q[:] = 0
P[:] = 0
q[idxmap[objfun.variable.value]] = 1
nothing
function processlinearterm(variables::Vector{VI}, coefficients::Vector, idxmap)
q = zeros(length(idxmap.varmap))
ncoeffs = length(coefficients)
@assert length(variables) == ncoeffs
for i = 1 : ncoeffs
q[idxmap[variables[i]].value] += coefficients[i]
end
q
end

function processconstraintfun!(A::AbstractMatrix, row::Int, idxmap, f::SingleVariable)
col = idxmap[f.variable].value
A[row, col] = 1
nothing
function symmetrize!(triplets::SparseTriplets)
I, J, V = triplets
n = length(V)
@assert length(I) == length(J) == n
for i = 1 : n
if I[i] != J[i]
push!(I, J[i])
push!(J, I[i])
push!(V, V[i])
end
end
end

function processconstraintfun!(A::AbstractMatrix, row::Int, idxmap, f::Affine)
n = length(f.coefficients)
@assert length(f.variables) == n
for i = 1 : n
var = f.variables[i]
coeff = f.coefficients[i]
col = idxmap[var].value
A[row, col] += coeff
function processconstraints(src::MOI.ModelLike, idxmap)
m = length(idxmap.conmap)
l = Vector{Float64}(uninitialized, m)
u = Vector{Float64}(uninitialized, m)
bounds = (l, u)

I = Int[]
J = Int[]
V = Float64[]
triplets = (I, J, V)

for (F, S) in MOI.get(src, MOI.ListOfConstraints())
processconstraints!(triplets, bounds, src, idxmap, F, S)
end
nothing
n = MOI.get(src, MOI.NumberOfVariables())
A = sparse(I, J, V, m, n)

(A, l, u)
end

function processconstraints!(A::SparseMatrixCSC, l::Vector, u::Vector, src::MOI.ModelLike, idxmap, F::Type, S::Type)
function processconstraints!(triplets::SparseTriplets, bounds::Tuple{<:Vector, <:Vector}, src::MOI.ModelLike, idxmap,
F::Type{<:MOI.AbstractFunction}, S::Type{<:MOI.AbstractSet})
cis_src = MOI.get(src, MOI.ListOfConstraintIndices{F, S}())
(l, u) = bounds
for ci in cis_src
s = MOI.Interval(MOI.get(src, MOI.ConstraintSet(), ci))
f = MOI.get(src, MOI.ConstraintFunction(), ci)
row = idxmap[ci].value
l[row] = s.lower - constant(f)
u[row] = s.upper - constant(f)
processconstraintfun!(A, row, idxmap, f)
processconstraints!(triplets, f, row, idxmap)
end
end

function processconstraints!(triplets::SparseTriplets, f::MOI.SingleVariable, row::Int, idxmap)
(I, J, V) = triplets
col = idxmap[f.variable].value
push!(I, row)
push!(J, col)
push!(V, 1)
nothing
end

function processconstraints!(triplets::SparseTriplets, f::MOI.ScalarAffineFunction, row::Int, idxmap)
(I, J, V) = triplets
ncoeff = length(f.coefficients)
@assert length(f.variables) == ncoeff
for i = 1 : ncoeff
var = f.variables[i]
coeff = f.coefficients[i]
col = idxmap[var].value
push!(I, row)
push!(J, col)
push!(V, coeff)
end
end

Expand All @@ -216,11 +242,12 @@ MOI.get(optimizer::OSQPOptimizer, ::MOI.ObjectiveSense) = optimizer.sense
MOI.canget(optimizer::OSQPOptimizer, ::MOI.NumberOfVariables) = !optimizer.isempty # https://github.com/oxfordcontrol/OSQP.jl/issues/10
MOI.get(optimizer::OSQPOptimizer, ::MOI.NumberOfVariables) = OSQP.dimensions(optimizer.model)[1]

MOI.canget(optimizer::OSQPOptimizer, ::MOI.ListOfVariableIndices) = MOI.get(optimizer, MOI.NumberOfVariables())
MOI.canget(optimizer::OSQPOptimizer, ::MOI.ListOfVariableIndices) = MOI.canget(optimizer, MOI.NumberOfVariables())
MOI.get(optimizer::OSQPOptimizer, ::MOI.ListOfVariableIndices) = [VI(i) for i = 1 : get(optimizer, MOI.NumberOfVariables())] # TODO: support for UnitRange would be nice

MOI.canget(optimizer::OSQPOptimizer, ::MOI.NumberOfConstraints) = !optimizer.isempty # https://github.com/oxfordcontrol/OSQP.jl/issues/10
MOI.get(optimizer::OSQPOptimizer, ::MOI.NumberOfConstraints) = OSQP.dimensions(optimizer.model)[2]
# FIXME or remove:
# MOI.canget(optimizer::OSQPOptimizer, ::MOI.NumberOfConstraints) = !optimizer.isempty # https://github.com/oxfordcontrol/OSQP.jl/issues/10
# MOI.get(optimizer::OSQPOptimizer, ::MOI.NumberOfConstraints) = OSQP.dimensions(optimizer.model)[2]

MOI.canget(optimizer::OSQPOptimizer, ::MOI.ListOfConstraints) = false # TODO
MOI.canget(optimizer::OSQPOptimizer, ::MOI.ListOfConstraintIndices) = false # TODO
Expand Down Expand Up @@ -380,7 +407,7 @@ end


## Constraints:
MOI.isvalid(optimizer::OSQPOptimizer, ci::CI) = ci.value 1 : get(optimizer, MOI.NumberOfConstraints())
# MOI.isvalid(optimizer::OSQPOptimizer, ci::CI) = ci.value ∈ 1 : get(optimizer, MOI.NumberOfConstraints()) # FIXME
MOI.canaddconstraint(optimizer::OSQPOptimizer, ::Type{F}, ::Type{S}) where {F <: MOI.AbstractFunction, S <: MOI.AbstractSet} = false

# TODO: can't modifyconstraint! with AbstractFunction because selective way to update A not exposed
Expand Down

0 comments on commit 44af411

Please sign in to comment.