From 0b1f98b8d2781eda18506cc25d0bb939a4aa10f3 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Tue, 9 Dec 2025 12:16:56 +0100 Subject: [PATCH 1/5] Simplify with MatrixOfConstraints --- src/matrix_input.jl | 95 ++++++++++++++++++++++----------------------- 1 file changed, 46 insertions(+), 49 deletions(-) diff --git a/src/matrix_input.jl b/src/matrix_input.jl index f0d5812..f5fc75d 100644 --- a/src/matrix_input.jl +++ b/src/matrix_input.jl @@ -79,66 +79,63 @@ function MOI.get( return _dot(model.A[ci.value, :]) end -""" - LPStandardForm{T, AT<:AbstractMatrix{T}, VT <: AbstractVector{T}} <: AbstractLPForm{T} +MOI.Utilities.@product_of_sets( + Equalities, + MOI.EqualTo{T}, +) -Represents a problem of the form: -``` -sense ⟨c, x⟩ -s.t. A x == b - x ≥ 0 -``` -""" -struct LPStandardForm{T,AT<:AbstractMatrix{T},VT<:AbstractVector{T}} <: - AbstractLPForm{T} - sense::MOI.OptimizationSense - c::VT - A::AT - b::VT -end -function MOI.get( - ::LPStandardForm{T}, - ::MOI.ListOfConstraintTypesPresent, -) where {T} - return [ - (MOI.ScalarAffineFunction{T}, MOI.EqualTo{T}), - (MOI.VectorOfVariables, MOI.Nonnegatives), - ] +# TODO specialize for SparseVector +function linear_function(c::AbstractVector{T}) where {T} + return MOI.ScalarAffineFunction( + [MOI.ScalarAffineTerm(c[i], MOI.VariableIndex(i)) for i in eachindex(c)], + zero(T), + ) end -const EQ{T} = MOI.ConstraintIndex{MOI.ScalarAffineFunction{T},MOI.EqualTo{T}} - -function MOI.get( - model::LPStandardForm{T}, - ::MOI.ListOfConstraintIndices{MOI.ScalarAffineFunction{T},MOI.EqualTo{T}}, -) where {T} - # TODO return `collect` with MOI v0.9.15 (see https://github.com/jump-dev/MathOptInterface.jl/pull/1110) - return collect(MOIU.LazyMap{EQ{T}}( - i -> EQ{T}(i), # FIXME `LazyMap` needs a `Function` so cannot just give `EQ{T}` - 1:size(model.A, 1), - )) +function linear_objective(sense::MOI.OptimizationSense, c::AbstractVector{T}) + model = MOI.Utilities.ObjectiveContainer{T}() + MOI.set(model, MOI.ObjectiveSense(), sense) + func = linear_function(c) + MOI.set(model, MOI.ObjectiveFunction{typeof(func)}(), func) + return model end -function MOI.get(model::LPStandardForm, ::MOI.ConstraintSet, ci::EQ) - return MOI.EqualTo(model.b[ci.value]) +function nonnegative_variables(n) + model = MOI.Utilities.VariablesContainer{T}() + for _ in 1:n + MOI.add_constrained_variable(model, MOI.GreaterThan(zero(T))) + end + return model end -const NONNEG = MOI.ConstraintIndex{MOI.VectorOfVariables,MOI.Nonnegatives} - -function MOI.get( - ::LPStandardForm, - ::MOI.ListOfConstraintIndices{MOI.VectorOfVariables,MOI.Nonnegatives}, -) - return [NONNEG(1)] +function equality_constraints(A::AbstractMatrix{T}, b::AbstractVector{T}) where {T} + return MOI.Utilities.MatrixOfConstraint{T}( + A, + b, + Equalities{T}(), # TODO + ) end -function MOI.get(model::LPStandardForm, ::MOI.ConstraintFunction, ci::NONNEG) - return MOI.VectorOfVariables(MOI.get(model, MOI.ListOfVariableIndices())) -end +""" + lp_standard_form(sense::MOI.OptimizationSense, c::AbtractVector, A::AbstractMatrix, b::AbstractVector) -function MOI.get(model::LPStandardForm, ::MOI.ConstraintSet, ci::NONNEG) - return MOI.Nonnegatives(MOI.get(model, MOI.NumberOfVariables())) +Represents a problem of the form: +``` +sense ⟨c, x⟩ +s.t. A x == b + x ≥ 0 +``` +""" +function lp_standard_form(sense::MOI.OptimizationSense, c::AbtractVector{T}, A::AbstractMatrix{T}, b::AbstractVector{T}) where {T} + m, n = size(A) + @assert length(c) == n + @assert length(b) == m + return MOI.Utilities.GenericModel{T}( + linear_objective(sense, c), + nonnegative_variables(n), + equality_constraints(A, b), + ) end """ From 7b697e245184ac625b6e7e79e56f93588211dd4d Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Tue, 9 Dec 2025 12:31:50 +0100 Subject: [PATCH 2/5] Fixes --- src/change_form.jl | 40 ---------------------------------------- src/matrix_input.jl | 6 +++--- 2 files changed, 3 insertions(+), 43 deletions(-) diff --git a/src/change_form.jl b/src/change_form.jl index 2b12b3d..c056437 100644 --- a/src/change_form.jl +++ b/src/change_form.jl @@ -15,21 +15,6 @@ function change_form(::Type{LPForm{T,AT,VT}}, lp::LPForm) where {T,AT,VT} ) end -function change_form( - ::Type{LPForm{T,AT,VT}}, - lp::LPStandardForm{T}, -) where {T,AT,VT} - return LPForm{T,AT,VT}( - lp.sense, - lp.c, - lp.A, - lp.b, - lp.b, - fill(zero(T), length(lp.c)), - fill(typemax(T), length(lp.c)), - ) -end - function change_form( ::Type{LPForm{T,AT,VT}}, lp::LPGeometricForm{T}, @@ -125,31 +110,6 @@ function change_form( return change_form(LPGeometricForm{T,AT,VT}, temp_lp) end -function change_form( - ::Type{LPStandardForm{T,AT,VT}}, - lp::LPStandardForm, -) where {T,AT,VT} - return LPStandardForm(lp.sense, lp.c, lp.A, lp.b) -end - -function change_form( - ::Type{LPStandardForm{T,AT,VT}}, - lp::LPGeometricForm{T}, -) where {T,AT,VT} - new_A = hcat(lp.A, -lp.A, AT(I, length(lp.b), length(lp.b))) - new_c = vcat(lp.c, -lp.c, fill(0.0, length(lp.b))) - return LPStandardForm{T,AT,VT}(lp.sense, new_c, new_A, copy(lp.b)) -end - -function change_form( - ::Type{LPStandardForm{T,AT,VT}}, - lp::F, -) where {T,AT,VT,F<:AbstractLPForm{T}} - temp_lp = change_form(LPForm{T,AT,VT}, lp) - new_lp = change_form(LPGeometricForm{T,AT,VT}, temp_lp) - return change_form(LPStandardForm{T,AT,VT}, new_lp) -end - function change_form( ::Type{LPSolverForm{T,AT,VT}}, lp::LPSolverForm, diff --git a/src/matrix_input.jl b/src/matrix_input.jl index f5fc75d..0c2ff21 100644 --- a/src/matrix_input.jl +++ b/src/matrix_input.jl @@ -93,7 +93,7 @@ function linear_function(c::AbstractVector{T}) where {T} ) end -function linear_objective(sense::MOI.OptimizationSense, c::AbstractVector{T}) +function linear_objective(sense::MOI.OptimizationSense, c::AbstractVector{T}) where {T} model = MOI.Utilities.ObjectiveContainer{T}() MOI.set(model, MOI.ObjectiveSense(), sense) func = linear_function(c) @@ -118,7 +118,7 @@ function equality_constraints(A::AbstractMatrix{T}, b::AbstractVector{T}) where end """ - lp_standard_form(sense::MOI.OptimizationSense, c::AbtractVector, A::AbstractMatrix, b::AbstractVector) + lp_standard_form(sense::MOI.OptimizationSense, c::AbstractVector, A::AbstractMatrix, b::AbstractVector) Represents a problem of the form: ``` @@ -127,7 +127,7 @@ s.t. A x == b x ≥ 0 ``` """ -function lp_standard_form(sense::MOI.OptimizationSense, c::AbtractVector{T}, A::AbstractMatrix{T}, b::AbstractVector{T}) where {T} +function lp_standard_form(sense::MOI.OptimizationSense, c::AbstractVector{T}, A::AbstractMatrix{T}, b::AbstractVector{T}) where {T} m, n = size(A) @assert length(c) == n @assert length(b) == m From ffb7cd54946a123c60cbd03787dab668bbc1000b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Fri, 19 Dec 2025 21:58:13 +0100 Subject: [PATCH 3/5] Fixes --- src/matrix_input.jl | 50 ++++++++++++++---- test/runtests.jl | 121 ++++++++++++++++++++++++++++---------------- 2 files changed, 119 insertions(+), 52 deletions(-) diff --git a/src/matrix_input.jl b/src/matrix_input.jl index 0c2ff21..8ffbf0e 100644 --- a/src/matrix_input.jl +++ b/src/matrix_input.jl @@ -84,6 +84,13 @@ MOI.Utilities.@product_of_sets( MOI.EqualTo{T}, ) +# We cannot use MOI.Utilities.Hyperrectangle in case the vector type is not `Vector` +struct Hyperrectangle{T,LT<:AbstractVector{T},UT<:AbstractVector{T}} <: MOI.Utilities.AbstractVectorBounds + lower::LT + upper::UT +end + +MOI.Utilities.function_constants(::Hyperrectangle{T}, row) where {T} = zero(T) # TODO specialize for SparseVector function linear_function(c::AbstractVector{T}) where {T} @@ -101,20 +108,44 @@ function linear_objective(sense::MOI.OptimizationSense, c::AbstractVector{T}) wh return model end -function nonnegative_variables(n) +function nonnegative_variables(n, ::Type{T}) where {T} model = MOI.Utilities.VariablesContainer{T}() for _ in 1:n - MOI.add_constrained_variable(model, MOI.GreaterThan(zero(T))) + vi = MOI.add_variable(model) + MOI.add_constraint(model, vi, MOI.GreaterThan(zero(T))) end return model end function equality_constraints(A::AbstractMatrix{T}, b::AbstractVector{T}) where {T} - return MOI.Utilities.MatrixOfConstraint{T}( - A, - b, - Equalities{T}(), # TODO - ) + sets = Equalities{T}() + for _ in eachindex(b) + MOI.Utilities.add_set(sets, MOI.Utilities.set_index(sets, MOI.EqualTo{T})) + end + constants = Hyperrectangle(b, b) + MOI.Utilities.final_touch(sets) + model = MOI.Utilities.MatrixOfConstraints{T}(A, constants, sets) + model.final_touch = true + return model +end + +# /!\ Type piracy +function MOI.Utilities.extract_function( + A::AbstractMatrix{T}, + row::Integer, + constant::T, +) where {T} + func = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm{T}[], constant) + for col in axes(A, 2) + val = A[row, col] + if !iszero(val) + push!( + func.terms, + MOI.ScalarAffineTerm(val, MOI.VariableIndex(col)), + ) + end + end + return func end """ @@ -131,11 +162,12 @@ function lp_standard_form(sense::MOI.OptimizationSense, c::AbstractVector{T}, A: m, n = size(A) @assert length(c) == n @assert length(b) == m - return MOI.Utilities.GenericModel{T}( + model = MOI.Utilities.GenericModel{T}( linear_objective(sense, c), - nonnegative_variables(n), + nonnegative_variables(n, T), equality_constraints(A, b), ) + return model end """ diff --git a/test/runtests.jl b/test/runtests.jl index 308ff0a..864071d 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -22,6 +22,7 @@ const dense_A = [ 1.0 2.0 3.0 4.0 ] + @testset "Matrix $(typeof(A))" for A in [dense_A, sparse(dense_A)] dense_b = [5.0, 6.0] dense_c = [7.0, 8.0] @@ -79,7 +80,7 @@ const dense_A = [ end @testset "change $(typeof(lp))" for lp in [ - MatOI.LPStandardForm{Float64,typeof(A),typeof(c)}( + MatOI.lp_standard_form( sense, c, A, @@ -105,13 +106,13 @@ const dense_A = [ ), ] test_expected(lp) - @testset "to $F" for F in [ - #MatOI.LPStandardForm{Float64, typeof(A)}, # FIXME doesn't work as the form gets bloated in the conversion - MatOI.LPForm{Float64,typeof(A),typeof(c)}, - MatOI.LPSolverForm{Float64,typeof(A),typeof(c)}, - ] - test_expected(MatOI.change_form(F, lp)) - end +# @testset "to $F" for F in [ +# MatOI.LPStandardForm{Float64, typeof(A)}, # FIXME doesn't work as the form gets bloated in the conversion +# MatOI.LPForm{Float64,typeof(A),typeof(c)}, +# MatOI.LPSolverForm{Float64,typeof(A),typeof(c)}, +# ] +# test_expected(MatOI.change_form(F, lp)) +# end end end @testset "Geometric form LP" begin @@ -157,41 +158,75 @@ const dense_A = [ ) end - @testset "change $(typeof(lp))" for lp in [ - MatOI.LPGeometricForm{Float64,typeof(A'),typeof(c)}( - sense, - b, - A', - c, - ), - MatOI.LPForm{Float64,typeof(A'),typeof(c)}( - sense, - b, - A', - [-Inf, -Inf], - c, - v_lb, - v_ub, - ), - MatOI.LPSolverForm{Float64,typeof(A'),typeof(c)}( - sense, - b, - A', - c, - [MatOI.LESS_THAN, MatOI.LESS_THAN], - v_lb, - v_ub, - ), - ] - test_expected(lp) - @testset "to $F" for F in [ - MatOI.LPGeometricForm{Float64,typeof(A),typeof(c)}, - MatOI.LPForm{Float64,typeof(A),typeof(c)}, - MatOI.LPSolverForm{Float64,typeof(A),typeof(c)}, - ] - test_expected(MatOI.change_form(F, lp)) - end - end +# @testset "change $(typeof(lp))" for lp in [ +# MatOI.LPGeometricForm{Float64,typeof(A'),typeof(c)}( +# sense, +# b, +# A', +# c, +# ), +# MatOI.LPForm{Float64,typeof(A'),typeof(c)}( +# sense, +# b, +# A', +# [-Inf, -Inf], +# c, +# v_lb, +# v_ub, +# ), +# MatOI.LPSolverForm{Float64,typeof(A'),typeof(c)}( +# sense, +# b, +# A', +# c, +# [MatOI.LESS_THAN, MatOI.LESS_THAN], +# v_lb, +# v_ub, +# ), +# ] +# test_expected(lp) +# @testset "to $F" for F in [ +# MatOI.LPGeometricForm{Float64,typeof(A),typeof(c)}, +# MatOI.LPForm{Float64,typeof(A),typeof(c)}, +# MatOI.LPSolverForm{Float64,typeof(A),typeof(c)}, +# ] +# test_expected(MatOI.change_form(F, lp)) +# end +# end end end end + +import MatrixOptInterface as MatOI +import MathOptInterface as MOI +m = 5 +n = 10 +c = rand(n) +b = rand(m) +A = rand(m, n) +model = MatOI.lp_standard_form(MOI.MIN_SENSE, c, A, b) +MOI.get(model, MOI.ListOfConstraintTypesPresent()) +println(model) + + +using JuMP, Clp +model = Model(Clp.Optimizer) +@variable(model, x) +@constraint(model, x == 1) +optimize!(model) +cache = backend(model).optimizer.model.model_cache.model +F = MOI.ScalarAffineFunction{Float64} +S = MOI.EqualTo{Float64} +ci = MOI.get(cache, MOI.ListOfConstraintIndices{F,S}())[] +mat = cache.constraints +r = MOI.Utilities.rows(mat, ci) +coef = mat.coefficients +cst = mat.constants +MOI.Utilities.extract_function( + coef, + r, + MOI.Utilities.function_constants(cst, r), +) +@edit MOI.Utilities.function_constants(cst, r) + +@edit MOI.get(mat, MOI.ConstraintFunction(), ci) From 2e698d573c763921a60103b60c0d0396957e99b6 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Fri, 19 Dec 2025 23:06:00 +0100 Subject: [PATCH 4/5] Clean up tests --- src/MatrixOptInterface.jl | 3 +- src/matrix_input.jl | 390 +++++++++++++++++++------------------- test/linear.jl | 189 ++++++++++++++++++ test/runtests.jl | 214 +-------------------- 4 files changed, 390 insertions(+), 406 deletions(-) create mode 100644 test/linear.jl diff --git a/src/MatrixOptInterface.jl b/src/MatrixOptInterface.jl index 0fbc923..951fde2 100644 --- a/src/MatrixOptInterface.jl +++ b/src/MatrixOptInterface.jl @@ -6,6 +6,7 @@ module MatrixOptInterface using LinearAlgebra, SparseArrays +import FillArrays using MathOptInterface const MOI = MathOptInterface const MOIU = MathOptInterface.Utilities @@ -67,6 +68,6 @@ end include("sparse_matrix.jl") include("conic_form.jl") include("matrix_input.jl") -include("change_form.jl") +#include("change_form.jl") end diff --git a/src/matrix_input.jl b/src/matrix_input.jl index 8ffbf0e..bf5209d 100644 --- a/src/matrix_input.jl +++ b/src/matrix_input.jl @@ -79,11 +79,6 @@ function MOI.get( return _dot(model.A[ci.value, :]) end -MOI.Utilities.@product_of_sets( - Equalities, - MOI.EqualTo{T}, -) - # We cannot use MOI.Utilities.Hyperrectangle in case the vector type is not `Vector` struct Hyperrectangle{T,LT<:AbstractVector{T},UT<:AbstractVector{T}} <: MOI.Utilities.AbstractVectorBounds lower::LT @@ -108,6 +103,67 @@ function linear_objective(sense::MOI.OptimizationSense, c::AbstractVector{T}) wh return model end +# /!\ Type piracy +function MOI.Utilities.extract_function( + A::Matrix{T}, + row::Integer, + constant::T, +) where {T} + func = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm{T}[], constant) + for col in axes(A, 2) + val = A[row, col] + if !iszero(val) + push!( + func.terms, + MOI.ScalarAffineTerm(val, MOI.VariableIndex(col)), + ) + end + end + return func +end + +# Copy-paste of MOI.Utilities.MutableSparseMatrixCSC +function _first_in_column( + A::SparseMatrixCSC, + row::Integer, + col::Integer, +) + range = SparseArrays.nzrange(A, col) + idx = searchsortedfirst(view(A.rowval, range), row) + return get(range, idx, last(range) + 1) +end + +# Copy-paste of MOI.Utilities.MutableSparseMatrixCSC +function MOI.Utilities.extract_function( + A::SparseMatrixCSC{T}, + row::Integer, + constant::T, +) where {T} + func = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm{T}[], constant) + for col in 1:A.n + idx = _first_in_column(A, row, col) + if idx > last(SparseArrays.nzrange(A, col)) + continue + end + r = A.rowval[idx] + if r == row + push!( + func.terms, + MOI.ScalarAffineTerm(A.nzval[idx], MOI.VariableIndex(col)), + ) + end + end + return func +end + +function free_variables(n, ::Type{T}) where {T} + model = MOI.Utilities.VariablesContainer{T}() + for _ in 1:n + MOI.add_variable(model) + end + return model +end + function nonnegative_variables(n, ::Type{T}) where {T} model = MOI.Utilities.VariablesContainer{T}() for _ in 1:n @@ -117,35 +173,87 @@ function nonnegative_variables(n, ::Type{T}) where {T} return model end +function interval_variables(lower::AbstractVector{T}, upper::AbstractVector{T}) where {T} + @assert eachindex(lower) == eachindex(upper) + model = MOI.Utilities.VariablesContainer{T}() + for i in eachindex(lower) + vi = MOI.add_variable(model) + MOI.add_constraint(model, vi, MOI.GreaterThan(lower[i])) + MOI.add_constraint(model, vi, MOI.LessThan(upper[i])) + end + return model +end + +MOI.Utilities.@product_of_sets( + EqualTos, + MOI.EqualTo{T}, +) + function equality_constraints(A::AbstractMatrix{T}, b::AbstractVector{T}) where {T} - sets = Equalities{T}() + sets = EqualTos{T}() for _ in eachindex(b) MOI.Utilities.add_set(sets, MOI.Utilities.set_index(sets, MOI.EqualTo{T})) end + MOI.Utilities.final_touch(sets) constants = Hyperrectangle(b, b) + model = MOI.Utilities.MatrixOfConstraints{T}(A, constants, sets) + model.final_touch = true + return model +end + +MOI.Utilities.@product_of_sets( + LessThans, + MOI.LessThan{T}, +) + +function lessthan_constraints(A::AbstractMatrix{T}, b::AbstractVector{T}) where {T} + sets = LessThans{T}() + for _ in eachindex(b) + MOI.Utilities.add_set(sets, MOI.Utilities.set_index(sets, MOI.LessThan{T})) + end MOI.Utilities.final_touch(sets) + constants = Hyperrectangle(FillArrays.Zeros{T}(length(b)), b) model = MOI.Utilities.MatrixOfConstraints{T}(A, constants, sets) model.final_touch = true return model end -# /!\ Type piracy -function MOI.Utilities.extract_function( - A::AbstractMatrix{T}, - row::Integer, - constant::T, -) where {T} - func = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm{T}[], constant) - for col in axes(A, 2) - val = A[row, col] - if !iszero(val) - push!( - func.terms, - MOI.ScalarAffineTerm(val, MOI.VariableIndex(col)), - ) - end +MOI.Utilities.@product_of_sets( + Intervals, + MOI.Interval{T}, +) + +function interval_constraints(A::AbstractMatrix{T}, lower::AbstractVector{T}, upper::AbstractVector{T}) where {T} + sets = Intervals{T}() + for _ in eachindex(lower) + MOI.Utilities.add_set(sets, MOI.Utilities.set_index(sets, MOI.Interval{T})) end - return func + MOI.Utilities.final_touch(sets) + constants = Hyperrectangle(lower, upper) + model = MOI.Utilities.MatrixOfConstraints{T}(A, constants, sets) + model.final_touch = true + return model +end + +MOI.Utilities.@mix_of_scalar_sets( + MixedLinearSets, + MOI.EqualTo{T}, + MOI.GreaterThan{T}, + MOI.LessThan{T}, +) + +function mix_of_constraints(A::AbstractMatrix{T}, b::AbstractVector{T}, senses::Vector{ConstraintSense}) where {T} + @assert eachindex(b) == eachindex(senses) + sets = MixedLinearSets{T}() + for sense in senses + # Int(EQUAL_TO) is 0 but `MOI.Utilities.set_index(sets, MOI.EqualTo{T})` is 1 + MOI.Utilities.add_set(sets, Int(sense) + 1) + end + MOI.Utilities.final_touch(sets) + constants = Hyperrectangle(b, b) + model = MOI.Utilities.MatrixOfConstraints{T}(A, constants, sets) + model.final_touch = true + return model end """ @@ -162,16 +270,15 @@ function lp_standard_form(sense::MOI.OptimizationSense, c::AbstractVector{T}, A: m, n = size(A) @assert length(c) == n @assert length(b) == m - model = MOI.Utilities.GenericModel{T}( + return MOI.Utilities.GenericModel{T}( linear_objective(sense, c), nonnegative_variables(n, T), equality_constraints(A, b), ) - return model end """ - LPGeometricForm{T, AT<:AbstractMatrix{T}, VT <: AbstractVector{T}} <: AbstractLPForm{T} + lp_geometric_form(sense::MOI.OptimizationSense, c::AbstractVector{T}, A::AbstractMatrix{T}, b::AbstractVector{T}) where {T} Represents a linear problem of the form: ``` @@ -179,133 +286,27 @@ sense ⟨c, x⟩ s.t. Ax <= b ``` """ -struct LPGeometricForm{T,AT<:AbstractMatrix{T},VT<:AbstractVector{T}} <: - AbstractLPForm{T} - sense::MOI.OptimizationSense - c::VT - A::AT - b::VT -end - -function MOI.get( - ::LPGeometricForm{T}, - ::MOI.ListOfConstraintTypesPresent, -) where {T} - return [(MOI.ScalarAffineFunction{T}, MOI.LessThan{T})] -end - -const LT{T} = MOI.ConstraintIndex{MOI.ScalarAffineFunction{T},MOI.LessThan{T}} - -function MOI.get( - model::LPGeometricForm{T}, - ::MOI.ListOfConstraintIndices{MOI.ScalarAffineFunction{T},MOI.LessThan{T}}, -) where {T} - # FIXME `copy_constraint` needs a `Vector` - return collect(MOIU.LazyMap{LT{T}}( - i -> LT{T}(i), # FIXME `LazyMap` needs a `Function` so cannot just give `LT{T}` - 1:size(model.A, 1), - )) -end - -function MOI.get(model::LPGeometricForm, ::MOI.ConstraintSet, ci::LT) - return MOI.LessThan(model.b[ci.value]) -end - -abstract type LPMixedForm{T} <: AbstractLPForm{T} end - -function MOI.get( - model::LPMixedForm{T}, - ::MOI.ListOfConstraintTypesPresent, -) where {T} - list = Tuple{DataType,DataType}[] - for S in - [MOI.EqualTo{T}, MOI.Interval{T}, MOI.GreaterThan{T}, MOI.LessThan{T}] - for F in [MOI.VariableIndex, MOI.ScalarAffineFunction{T}] - if !iszero(MOI.get(model, MOI.NumberOfConstraints{F,S}())) - push!(list, (F, S)) - end - end - end - return list -end - -struct BoundSense <: MOI.AbstractVariableAttribute end - -function MOI.get(model::LPMixedForm, ::BoundSense, vi::MOI.VariableIndex) - return _bound_sense(model.v_lb[vi.value], model.v_ub[vi.value]) -end - -const LinearBounds{T} = - Union{MOI.EqualTo{T},MOI.Interval{T},MOI.GreaterThan{T},MOI.LessThan{T}} - -function MOI.get( - model::LPMixedForm{T}, - ::MOI.NumberOfConstraints{MOI.ScalarAffineFunction{T},S}, -) where {T,S<:LinearBounds{T}} - s = _sense(S) - return count(1:size(model.A, 1)) do i - return _constraint_bound_sense(model, i) == s - end -end - -const AFF{T,S} = MOI.ConstraintIndex{MOI.ScalarAffineFunction{T},S} - -function MOI.get( - model::LPMixedForm{T}, - ::MOI.ListOfConstraintIndices{MOI.ScalarAffineFunction{T},S}, -) where {T,S<:LinearBounds{T}} - s = _sense(S) - # FIXME `copy_constraint` needs a `Vector` - return collect( - MOIU.LazyMap{AFF{T,S}}( - i -> AFF{T,S}(i), # FIXME `LazyMap` needs a `Function` so cannot just give `AFF{T, S}` - Base.Iterators.Filter(1:size(model.A, 1)) do i - return _constraint_bound_sense(model, i) == s - end, - ), - ) -end - -const VBOUND{S} = MOI.ConstraintIndex{MOI.VariableIndex,S} - -function MOI.get( - model::LPMixedForm{T}, - ::MOI.NumberOfConstraints{MOI.VariableIndex,S}, -) where {T,S<:LinearBounds{T}} - s = _sense(S) - return count(MOI.get(model, MOI.ListOfVariableIndices())) do vi - return MOI.get(model, BoundSense(), vi) == s - end -end - -function MOI.get( - model::LPMixedForm{T}, - ::MOI.ListOfConstraintIndices{MOI.VariableIndex,S}, -) where {T,S<:LinearBounds{T}} - s = _sense(S) - return collect( - MOIU.LazyMap{VBOUND{S}}( - Base.Iterators.Filter( - MOI.get(model, MOI.ListOfVariableIndices()), - ) do vi - return MOI.get(model, BoundSense(), vi) == s - end, - ) do vi - return VBOUND{S}(vi.value) - end, +function lp_geometric_form(sense::MOI.OptimizationSense, c::AbstractVector{T}, A::AbstractMatrix{T}, b::AbstractVector{T}) where {T} + m, n = size(A) + @assert length(c) == n + @assert length(b) == m + return MOI.Utilities.GenericModel{T}( + linear_objective(sense, c), + free_variables(n, T), + lessthan_constraints(A, b), ) end -function MOI.get(::LPMixedForm, ::MOI.ConstraintFunction, ci::VBOUND) - return MOI.VariableIndex(ci.value) -end - -function MOI.get(model::LPMixedForm, ::MOI.ConstraintSet, ci::VBOUND) - return _bound_set(model.v_lb[ci.value], model.v_ub[ci.value]) -end - """ - LPForm{T, AT<:AbstractMatrix{T}, VT <: AbstractVector{T}} + function lp_form( + sense::MOI.OptimizationSense, + c::AbstractVector{T}, + A::AbstractMatrix{T}, + c_lb::AbstractVector{T}, + c_ub::AbstractVector{T}, + v_lb::AbstractVector{T}, + v_ub::AbstractVector{T}, + ) where {T} Represents a problem of the form: ``` @@ -314,26 +315,38 @@ s.t. c_lb <= Ax <= c_ub v_lb <= x <= v_ub ``` """ -struct LPForm{T,AT<:AbstractMatrix{T},VT<:AbstractVector{T}} <: LPMixedForm{T} #, V<:AbstractVector{T} #, M<:AbstractMatrix{T}} - sense::MOI.OptimizationSense - c::VT - A::AT - c_lb::VT - c_ub::VT - v_lb::VT - v_ub::VT -end - -function _constraint_bound_sense(model::LPForm, i) - return _bound_sense(model.c_lb[i], model.c_ub[i]) -end - -function MOI.get(model::LPForm, ::MOI.ConstraintSet, ci::AFF) - return _bound_set(model.c_lb[ci.value], model.c_ub[ci.value]) +function lp_form( + sense::MOI.OptimizationSense, + c::AbstractVector{T}, + A::AbstractMatrix{T}, + c_lb::AbstractVector{T}, + c_ub::AbstractVector{T}, + v_lb::AbstractVector{T}, + v_ub::AbstractVector{T}, +) where {T} + m, n = size(A) + @assert length(c) == n + @assert length(v_lb) == n + @assert length(v_ub) == n + @assert length(c_lb) == m + @assert length(c_ub) == m + return MOI.Utilities.GenericModel{T}( + linear_objective(sense, c), + interval_variables(v_lb, v_ub), + interval_constraints(A, c_lb, c_ub), + ) end """ - LPSolverForm{T, AT<:AbstractMatrix{T}, VT<:AbstractVector{T}} + lp_solver_form( + sense::MOI.OptimizationSense, + c::AbstractVector{T}, + A::AbstractMatrix{T}, + b::AbstractVector{T}, + sense::Vector{ConstraintSense}, + v_lb::AbstractVector{T}, + v_ub::AbstractVector{T}, + ) where {T} Represents a problem of the form: ``` @@ -342,32 +355,25 @@ s.t. Ax senses b v_lb <= x <= v_ub ``` """ -struct LPSolverForm{T,AT<:AbstractMatrix{T},VT<:AbstractVector{T}} <: - LPMixedForm{T} - sense::MOI.OptimizationSense - c::VT - A::AT - b::VT - senses::Vector{ConstraintSense} - v_lb::VT - v_ub::VT -end - -function _constraint_bound_sense(model::LPSolverForm, i) - return model.senses[i] -end - -function MOI.get(model::LPSolverForm, ::MOI.ConstraintSet, ci::AFF) - s = model.senses[ci.value] - β = model.b[ci.value] - if s == GREATER_THAN - return MOI.GreaterThan(β) - elseif s == LESS_THAN - return MOI.LessThan(β) - else - s == EQUAL_TO || error("Invalid $s") - return MOI.EqualTo(β) - end +function lp_solver_form( + sense::MOI.OptimizationSense, + c::AbstractVector{T}, + A::AbstractMatrix{T}, + b::AbstractVector{T}, + senses::Vector{ConstraintSense}, + v_lb::AbstractVector{T}, + v_ub::AbstractVector{T}, +) where {T} + m, n = size(A) + @assert length(c) == n + @assert length(v_lb) == n + @assert length(v_ub) == n + @assert length(b) == m + return MOI.Utilities.GenericModel{T}( + linear_objective(sense, c), + interval_variables(v_lb, v_ub), + mix_of_constraints(A, b, senses), + ) end """ @@ -376,7 +382,7 @@ end A mixed-integer problem represented by a linear problem of type `LP` and a vector indicating each `VariableType`. """ -struct MILP{T,LP<:AbstractLPForm{T}} +struct MILP{T,LP<:MOI.ModelLike} lp::LP variable_type::Vector{VariableType} end diff --git a/test/linear.jl b/test/linear.jl new file mode 100644 index 0000000..3eda2c1 --- /dev/null +++ b/test/linear.jl @@ -0,0 +1,189 @@ +module TestLinear + +using SparseArrays, Test +import MathOptInterface as MOI +import MatrixOptInterface as MatOI + +function runtests() + for name in names(@__MODULE__; all = true) + if startswith("$(name)", "test_") + @testset "$(name)" begin + getfield(@__MODULE__, name)() + end + end + end + return +end + +function _test_expected(form, expected, var_names, con_names, S) + model = MOI.Utilities.Model{Float64}() + MOI.copy_to( + MOI.Bridges.Constraint.Scalarize{Float64}(model), + form, + ) + MOI.set( + model, + MOI.VariableName(), + MOI.VariableIndex.(eachindex(var_names)), + var_names, + ) + MOI.set( + model, + MOI.ConstraintName(), + MOI.ConstraintIndex{ + MOI.ScalarAffineFunction{Float64}, + S, + }.(eachindex(con_names)), + con_names, + ) + return MOI.Test.util_test_models_equal( + model, + expected, + var_names, + con_names, + ) +end + + +function test_standard_form() + dense_A = [ + 1.0 2.0 + 3.0 4.0 + ] + dense_b = [5.0, 6.0] + dense_c = [7.0, 8.0] + s = """ + variables: x1, x2 + cx1: x1 >= 0.0 + cx2: x2 >= 0.0 + c1: x1 + 2x2 == 5.0 + c2: 3x1 + 4x2 == 6.0 + minobjective: 7x1 + 8x2 + """ + expected = MOI.Utilities.Model{Float64}() + MOI.Utilities.loadfromstring!(expected, s) + + var_names = ["x1", "x2"] + con_names = ["c1", "c2"] + + sense = MOI.MIN_SENSE + v_lb = [0.0, 0.0] + v_ub = [Inf, Inf] + + @testset "Matrix $(typeof(A))" for A in [dense_A, sparse(dense_A)] + @testset "Vector $(typeof(b))" for (b, c) in zip( + [dense_b, sparsevec(dense_b)], + [dense_c, sparsevec(dense_c)], + ) + @testset "change $func" for (func, args) in [ + (MatOI.lp_standard_form, ( + sense, + c, + A, + b, + )), + (MatOI.lp_solver_form, ( + sense, + c, + A, + b, + [MatOI.EQUAL_TO, MatOI.EQUAL_TO], + v_lb, + v_ub, + )), + ] + lp = func(args...) + _test_expected(lp, expected, var_names, con_names, MOI.EqualTo{Float64}) + end + end + end +end + +function test_geometric_form() + dense_A = [ + 1.0 2.0 + 3.0 4.0 + ] + dense_b = [5.0, 6.0] + dense_c = [7.0, 8.0] + s = """ + variables: x1, x2 + cx1: x1 >= 0.0 + cx2: x2 >= 0.0 + c1: x1 + 2x2 <= 5.0 + c2: 3x1 + 4x2 <= 6.0 + minobjective: 7x1 + 8x2 + """ + expected = MOI.Utilities.Model{Float64}() + MOI.Utilities.loadfromstring!(expected, s) + + var_names = ["x1", "x2"] + con_names = ["c1", "c2"] + + sense = MOI.MIN_SENSE + + @testset "Matrix $(typeof(A))" for A in [dense_A, sparse(dense_A)] + @testset "Vector $(typeof(b))" for (b, c) in zip( + [dense_b, sparsevec(dense_b)], + [dense_c, sparsevec(dense_c)], + ) + lp = MatOI.lp_geometric_form( + sense, + c, + A, + b, + ) + _test_expected(lp, expected, var_names, con_names, MOI.LessThan{Float64}) + end + end +end + +function test_form() + dense_A = [ + 1.0 2.0 + 3.0 4.0 + ] + dense_lb = [5.0, 6.0] + dense_ub = [7.0, 9.0] + dense_c = [7.0, 8.0] + s = """ + variables: x1, x2 + cx1: x1 >= 0.0 + cx2: x2 >= 0.0 + c1: x1 + 2x2 in MOI.Interval(5.0, 7.0) + c2: 3x1 + 4x2 in MOI.Interval(6.0, 9.0) + minobjective: 7x1 + 8x2 + """ + expected = MOI.Utilities.Model{Float64}() + MOI.Utilities.loadfromstring!(expected, s) + + var_names = ["x1", "x2"] + con_names = ["c1", "c2"] + + sense = MOI.MIN_SENSE + v_lb = [0.0, 0.0] + v_ub = [Inf, Inf] + + @testset "Matrix $(typeof(A))" for A in [dense_A, sparse(dense_A)] + @testset "Vector $(typeof(c))" for (c_lb, c_ub, c) in zip( + [dense_lb, sparsevec(dense_lb)], + [dense_ub, sparsevec(dense_ub)], + [dense_c, sparsevec(dense_c)], + ) + lp = MatOI.lp_form( + sense, + c, + A, + c_lb, + c_ub, + v_lb, + v_ub, + ) + _test_expected(lp, expected, var_names, con_names, MOI.Interval{Float64}) + end + end +end + +end + +TestLinear.runtests() diff --git a/test/runtests.jl b/test/runtests.jl index 864071d..d95aa47 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -16,217 +16,5 @@ const MOIT = MOI.Test const ATOL = 1e-4 const RTOL = 1e-4 +include("linear.jl") include("conic_form.jl") - -const dense_A = [ - 1.0 2.0 - 3.0 4.0 -] - -@testset "Matrix $(typeof(A))" for A in [dense_A, sparse(dense_A)] - dense_b = [5.0, 6.0] - dense_c = [7.0, 8.0] - @testset "Vector $(typeof(b))" for (b, c) in zip( - [dense_b, sparsevec(dense_b)], - [dense_c, sparsevec(dense_c)], - ) - @testset "Standard form LP" begin - s = """ - variables: x1, x2 - cx1: x1 >= 0.0 - cx2: x2 >= 0.0 - c1: x1 + 2x2 == 5.0 - c2: 3x1 + 4x2 == 6.0 - minobjective: 7x1 + 8x2 - """ - expected = MOIU.Model{Float64}() - MOIU.loadfromstring!(expected, s) - - var_names = ["x1", "x2"] - con_names = ["c1", "c2"] - vcon_names = ["cx1", "cx2"] - model = MOIU.Model{Float64}() - - sense = MOI.MIN_SENSE - v_lb = [0.0, 0.0] - v_ub = [Inf, Inf] - - function test_expected(form) - MOI.copy_to( - MOI.Bridges.Constraint.Scalarize{Float64}(model), - form, - ) - MOI.set( - model, - MOI.VariableName(), - MOI.VariableIndex.(1:2), - var_names, - ) - MOI.set( - model, - MOI.ConstraintName(), - MOI.ConstraintIndex{ - MOI.ScalarAffineFunction{Float64}, - MOI.EqualTo{Float64}, - }.(1:2), - con_names, - ) - return MOIT.util_test_models_equal( - model, - expected, - var_names, - con_names, - ) - end - - @testset "change $(typeof(lp))" for lp in [ - MatOI.lp_standard_form( - sense, - c, - A, - b, - ), - MatOI.LPForm{Float64,typeof(A),typeof(c)}( - sense, - c, - A, - b, - b, - v_lb, - v_ub, - ), - MatOI.LPSolverForm{Float64,typeof(A),typeof(c)}( - sense, - c, - A, - b, - [MatOI.EQUAL_TO, MatOI.EQUAL_TO], - v_lb, - v_ub, - ), - ] - test_expected(lp) -# @testset "to $F" for F in [ -# MatOI.LPStandardForm{Float64, typeof(A)}, # FIXME doesn't work as the form gets bloated in the conversion -# MatOI.LPForm{Float64,typeof(A),typeof(c)}, -# MatOI.LPSolverForm{Float64,typeof(A),typeof(c)}, -# ] -# test_expected(MatOI.change_form(F, lp)) -# end - end - end - @testset "Geometric form LP" begin - s = """ - variables: x1, x2 - c1: x1 + 3x2 <= 7.0 - c2: 2x1 + 4x2 <= 8.0 - maxobjective: 5x1 + 6x2 - """ - expected = MOIU.Model{Float64}() - MOIU.loadfromstring!(expected, s) - - var_names = ["x1", "x2"] - con_names = ["c1", "c2"] - model = MOIU.Model{Float64}() - - sense = MOI.MAX_SENSE - v_lb = [-Inf, -Inf] - v_ub = [Inf, Inf] - - function test_expected(form) - MOI.copy_to(model, form) - MOI.set( - model, - MOI.VariableName(), - MOI.VariableIndex.(1:2), - var_names, - ) - MOI.set( - model, - MOI.ConstraintName(), - MOI.ConstraintIndex{ - MOI.ScalarAffineFunction{Float64}, - MOI.LessThan{Float64}, - }.(1:2), - con_names, - ) - return MOIT.util_test_models_equal( - model, - expected, - var_names, - con_names, - ) - end - -# @testset "change $(typeof(lp))" for lp in [ -# MatOI.LPGeometricForm{Float64,typeof(A'),typeof(c)}( -# sense, -# b, -# A', -# c, -# ), -# MatOI.LPForm{Float64,typeof(A'),typeof(c)}( -# sense, -# b, -# A', -# [-Inf, -Inf], -# c, -# v_lb, -# v_ub, -# ), -# MatOI.LPSolverForm{Float64,typeof(A'),typeof(c)}( -# sense, -# b, -# A', -# c, -# [MatOI.LESS_THAN, MatOI.LESS_THAN], -# v_lb, -# v_ub, -# ), -# ] -# test_expected(lp) -# @testset "to $F" for F in [ -# MatOI.LPGeometricForm{Float64,typeof(A),typeof(c)}, -# MatOI.LPForm{Float64,typeof(A),typeof(c)}, -# MatOI.LPSolverForm{Float64,typeof(A),typeof(c)}, -# ] -# test_expected(MatOI.change_form(F, lp)) -# end -# end - end - end -end - -import MatrixOptInterface as MatOI -import MathOptInterface as MOI -m = 5 -n = 10 -c = rand(n) -b = rand(m) -A = rand(m, n) -model = MatOI.lp_standard_form(MOI.MIN_SENSE, c, A, b) -MOI.get(model, MOI.ListOfConstraintTypesPresent()) -println(model) - - -using JuMP, Clp -model = Model(Clp.Optimizer) -@variable(model, x) -@constraint(model, x == 1) -optimize!(model) -cache = backend(model).optimizer.model.model_cache.model -F = MOI.ScalarAffineFunction{Float64} -S = MOI.EqualTo{Float64} -ci = MOI.get(cache, MOI.ListOfConstraintIndices{F,S}())[] -mat = cache.constraints -r = MOI.Utilities.rows(mat, ci) -coef = mat.coefficients -cst = mat.constants -MOI.Utilities.extract_function( - coef, - r, - MOI.Utilities.function_constants(cst, r), -) -@edit MOI.Utilities.function_constants(cst, r) - -@edit MOI.get(mat, MOI.ConstraintFunction(), ci) From 062cc1f9fafd7c5b7604e1dbf6d57070cde5e0cf Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Fri, 19 Dec 2025 23:06:34 +0100 Subject: [PATCH 5/5] Add FillArrays --- Project.toml | 2 ++ 1 file changed, 2 insertions(+) diff --git a/Project.toml b/Project.toml index fa7b579..5323233 100644 --- a/Project.toml +++ b/Project.toml @@ -3,11 +3,13 @@ uuid = "2f4eb8e6-3e35-4ae4-8c7a-f3d7d9bf20ed" version = "0.1.0" [deps] +FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" [compat] +FillArrays = "1.15.0" MathOptInterface = "1" julia = "1.6"