diff --git a/Project.toml b/Project.toml index 3ae697e..9568ccb 100644 --- a/Project.toml +++ b/Project.toml @@ -12,6 +12,7 @@ SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" [compat] Graphs = "1" JuMP = "1" +MathOptInterface = "1" julia = "1" [extras] diff --git a/src/identify_variables.jl b/src/identify_variables.jl index e757553..f28ed62 100644 --- a/src/identify_variables.jl +++ b/src/identify_variables.jl @@ -28,12 +28,6 @@ import MathOptInterface as MOI import MathProgIncidence: get_equality_constraints -# TODO: This file implements functions that filter duplicates from the -# vectors of identified variables. It may be useful at somepoint to -# identify variables, preserving duplicates. This can be implemented -# when/if there is a need. - - """ identify_unique_variables(constraints::Vector)::Vector{JuMP.VariableRef} @@ -207,50 +201,149 @@ function identify_unique_variables( return _filter_duplicates(refs) end - """ identify_unique_variables(fcn)::Vector{JuMP.VariableIndex} Return the variables that appear in the provided MathOptInterface function. +No variable will appear more than once. + # Implementation -Only `ScalarQuadraticFunction` and `ScalarAffineFunction` are supported. -This can be changed there is demand for other functions. For each type of -supported function, the `_get_variable_terms` function should be defined. -Then, for the type of each term, an additional `identify_unique_variables` -function should be implemented. +Only `ScalarNonlinearFunction`, `ScalarQuadraticFunction`, and +`ScalarAffineFunction` are supported. This can be changed there is demand for +other functions. +These methods are implemented by first identifying all variables that +participate in the function, then filtering out duplicate variables. """ function identify_unique_variables( - fcn::Union{MOI.ScalarQuadraticFunction, MOI.ScalarAffineFunction}, + fcn::Union{ + MOI.ScalarAffineFunction, + MOI.ScalarQuadraticFunction, + MOI.ScalarNonlinearFunction, + }, )::Vector{MOI.VariableIndex} - variables = Vector{MOI.VariableIndex}() - for terms in _get_variable_terms(fcn) - for term in terms - for var in identify_unique_variables(term) - push!(variables, var) - end - end - end + variables = _identify_variables(fcn) return _filter_duplicates(variables) end +# This method is used to handle function-in-set constraints. +function identify_unique_variables( + var::MOI.VariableIndex +)::Vector{MOI.VariableIndex} + return [var] +end + +# NOTE: We will get a MethodError if this is called with a non-vector, +# non-Affine/Quadratic/Nonlinear function. Should probably implement +# some default method to catch that. function identify_unique_variables( - fcn::T -)::Vector{MOI.VariableIndex} where {T<:MOI.AbstractVectorFunction} + fcn::MOI.AbstractVectorFunction +)::Vector{MOI.VariableIndex} throw(TypeError( fcn, - Union{MOI.ScalarQuadraticFunction, MOI.ScalarAffineFunction}, + Union{ + MOI.ScalarAffineFunction, + MOI.ScalarQuadraticFunction, + MOI.ScalarNonlinearFunction, + }, typeof(fcn), )) end -function identify_unique_variables( - var::MOI.VariableIndex +""" + _identify_variables(fcn)::Vector{JuMP.VariableIndex} + +Return all variables that appear in the provided MathOptInterface function. + +Duplicates may be present in the resulting vector of variables. +If there is a use case, these methods can be made public. + +# Implementation +Implemented for `ScalarAffineFunction`, `ScalarQuadraticFunction`, and +`ScalarNonlinearFunction`. Relies on an underlying _collect_variables! +method that (potentially recursively) builds up a vector of variables +in-place. + +""" +function _identify_variables( + fcn::Union{ + MOI.ScalarAffineFunction, + MOI.ScalarQuadraticFunction, + MOI.ScalarNonlinearFunction, + }, )::Vector{MOI.VariableIndex} - return [var] + variables = Vector{MOI.VariableIndex}() + _collect_variables!(variables, fcn) + return variables end +""" + _collect_variables!(variables, fcn)::Vector{JuMP.VariableIndex} + +Add variables from `fcn` to the `variables` vector + +# Implementation +Implemented for `ScalarAffineFunction`, `ScalarQuadraticFunction`, and +`ScalarNonlinearFunction`. For affine and quadratic functions, we iterate +over terms with the `_get_variable_terms` function. For nonlinear functions, +we recurse until pushing a root node onto the variable stack (or hitting +an affine/quadratic function). Methods may need to be added if more node types +are added to `ScalarNonlinearFunction` in the future. + +""" +function _collect_variables!( + variables::Vector{MOI.VariableIndex}, + fcn::MOI.ScalarNonlinearFunction, +)::Vector{MOI.VariableIndex} + for arg in fcn.args + _collect_variables!(variables, arg) + end + return variables +end + +function _collect_variables!( + variables::Vector{MOI.VariableIndex}, + fcn::Union{MOI.ScalarQuadraticFunction, MOI.ScalarAffineFunction}, +)::Vector{MOI.VariableIndex} + for terms in _get_variable_terms(fcn) + for term in terms + _collect_variables!(variables, term) + end + end + return variables +end + +function _collect_variables!( + variables::Vector{MOI.VariableIndex}, + var::MOI.VariableIndex, +)::Vector{MOI.VariableIndex} + push!(variables, var) + return variables +end + +function _collect_variables!( + variables::Vector{MOI.VariableIndex}, + var::Float64, +)::Vector{MOI.VariableIndex} + return variables +end + +function _collect_variables!( + variables::Vector{MOI.VariableIndex}, + term::MOI.ScalarAffineTerm, +)::Vector{MOI.VariableIndex} + push!(variables, term.variable) + return variables +end + +function _collect_variables!( + variables::Vector{MOI.VariableIndex}, + term::MOI.ScalarQuadraticTerm, +)::Vector{MOI.VariableIndex} + push!(variables, term.variable_1, term.variable_2) + return variables +end """ _get_variable_terms(fcn) @@ -301,7 +394,6 @@ function identify_unique_variables( end end - """ _filter_duplicates @@ -324,7 +416,6 @@ function _filter_duplicates( return filtered end - function _filter_duplicates( variables::Vector{JuMP.VariableRef}, )::Vector{JuMP.VariableRef} diff --git a/test/identify_variables.jl b/test/identify_variables.jl index a132c21..42ea01f 100644 --- a/test/identify_variables.jl +++ b/test/identify_variables.jl @@ -26,8 +26,8 @@ using MathProgIncidence: identify_unique_variables # Local import of JuMP models for testing include("models.jl") # make_degenerate_flow_model, make_simple_model -function test_linear() - m = make_degenerate_flow_model() +function test_linear(model_function=make_degenerate_flow_model) + m = model_function() variables = identify_unique_variables(m[:sum_comp_eqn]) # What is happening when we hash a VariableRef? I.e. how does # the model get hashed? @@ -36,16 +36,16 @@ function test_linear() @test(pred_var_set == Set(variables)) end -function test_quadratic() - m = make_degenerate_flow_model() +function test_quadratic(model_function=make_degenerate_flow_model) + m = model_function() variables = identify_unique_variables(m[:comp_dens_eqn][1]) pred_var_set = Set([m[:rho], m[:x][1]]) @test(length(variables) == length(pred_var_set)) @test(pred_var_set == Set(variables)) end -function test_nonlinear() - m = make_degenerate_flow_model() +function test_nonlinear(model_function=make_degenerate_flow_model) + m = model_function() variables = identify_unique_variables(m[:bulk_dens_eqn]) pred_var_set = Set([m[:rho], m[:x][1], m[:x][2], m[:x][3]]) @test(length(variables) == length(pred_var_set)) @@ -61,8 +61,8 @@ function test_nonlinear_with_potential_duplicates() @test(Set(variables) == Set([m[:var][1], m[:var][2]])) end -function test_several_equalities() - m = make_degenerate_flow_model() +function test_several_equalities(model_function=make_degenerate_flow_model) + m = model_function() constraints = [m[:comp_flow_eqn][1], m[:sum_comp_eqn], m[:bulk_dens_eqn]] variables = identify_unique_variables(constraints) pred_var_set = Set([ @@ -72,8 +72,8 @@ function test_several_equalities() @test(pred_var_set == Set(variables)) end -function test_several_constraints_with_ineq() - m = make_degenerate_flow_model() +function test_several_constraints_with_ineq(model_function=make_degenerate_flow_model) + m = model_function() @JuMP.constraint(m, ineq1, m[:flow_comp][2] >= 1.0) @JuMP.constraint(m, ineq2, m[:flow_comp][1]^2 + m[:flow_comp][3]^2 <= 1.0) constraints = [ @@ -98,8 +98,8 @@ function test_several_constraints_with_ineq() @test(pred_var_set == Set(variables)) end -function test_model() - m = make_degenerate_flow_model() +function test_model(model_function=make_degenerate_flow_model) + m = model_function() @JuMP.variable(m, dummy) @JuMP.NLconstraint(m, dummy^3.0 <= 5) # Note that include_inequalities=false by default. @@ -118,8 +118,8 @@ function test_model() @test(pred_var_set == Set(variables)) end -function test_model_with_ineq() - m = make_degenerate_flow_model() +function test_model_with_ineq(model_function=make_degenerate_flow_model) + m = model_function() @JuMP.variable(m, dummy) @JuMP.NLconstraint(m, dummy^3.0 <= 5) variables = identify_unique_variables(m, include_inequality=true) @@ -149,8 +149,8 @@ function test_function_with_variable_squared() @test(Set(variables) == Set([m[:dummy1].index, m[:dummy2].index])) end -function test_model_bad_constr() - m = make_degenerate_flow_model() +function test_model_bad_constr(model_function=make_degenerate_flow_model) + m = model_function() @JuMP.variable(m, dummy[1:2]) @JuMP.constraint(m, vectorcon, dummy in MOI.Nonnegatives(2)) @test_throws( @@ -159,8 +159,8 @@ function test_model_bad_constr() ) end -function test_model_bad_constr_no_ineq() - m = make_degenerate_flow_model() +function test_model_bad_constr_no_ineq(model_function=make_degenerate_flow_model) + m = model_function() @JuMP.variable(m, dummy[1:2]) @JuMP.constraint(m, vectorcon, dummy in MOI.Nonnegatives(2)) # Note that we don't throw an error because we don't attempt to @@ -208,8 +208,8 @@ function test_two_constraints_same_type() @test Set(vars) == pred_var_set end -function test_variables_in_inequalities() - m = make_simple_model() +function test_variables_in_inequalities(model_function=make_simple_model) + m = model_function() @JuMP.variable(m, y >= 1) x = m[:x] @JuMP.objective(m, Min, x[1] + 2*x[2] + 3*x[3] + y^2) @@ -224,18 +224,39 @@ end @testset "get-equality" begin test_linear() + test_linear(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_quadratic() + test_quadratic(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_nonlinear() + test_nonlinear(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_nonlinear_with_potential_duplicates() + test_several_equalities() + test_several_equalities(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_several_constraints_with_ineq() + test_several_constraints_with_ineq(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_model() + test_model(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_model_with_ineq() + test_model_with_ineq(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_model_bad_constr() + test_model_bad_constr(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_model_bad_constr_no_ineq() + test_model_bad_constr_no_ineq(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_function_with_variable_squared() test_fixing_constraint() test_inequality_with_bounds() test_two_constraints_same_type() + test_variables_in_inequalities() + test_variables_in_inequalities(make_simple_model_with_ScalarNonlinearFunction) end diff --git a/test/interface.jl b/test/interface.jl index e7aed69..6af66c2 100644 --- a/test/interface.jl +++ b/test/interface.jl @@ -41,8 +41,8 @@ function _test_igraph_fields(igraph, constraints, variables) end -function test_construct_interface() - m = make_degenerate_flow_model() +function test_construct_interface(model_function=make_degenerate_flow_model) + m = model_function() igraph = MathProgIncidence.IncidenceGraphInterface(m) variables = [ @@ -70,8 +70,8 @@ function test_construct_interface() end -function test_construct_interface_rectangular() - m = make_degenerate_flow_model() +function test_construct_interface_rectangular(model_function=make_degenerate_flow_model) + m = model_function() @JuMP.constraint( m, sum_flow_eqn, @@ -105,8 +105,8 @@ function test_construct_interface_rectangular() end -function test_get_adjacent_to_linear_constraint() - m = make_degenerate_flow_model() +function test_get_adjacent_to_linear_constraint(model_function=make_degenerate_flow_model) + m = model_function() igraph = MathProgIncidence.IncidenceGraphInterface(m) con = m[:sum_comp_eqn] adjacent = MathProgIncidence.get_adjacent(igraph, con) @@ -115,8 +115,8 @@ function test_get_adjacent_to_linear_constraint() end -function test_get_adjacent_to_quadratic_constraint() - m = make_degenerate_flow_model() +function test_get_adjacent_to_quadratic_constraint(model_function=make_degenerate_flow_model) + m = model_function() igraph = MathProgIncidence.IncidenceGraphInterface(m) con = m[:comp_dens_eqn][1] adjacent = MathProgIncidence.get_adjacent(igraph, con) @@ -125,8 +125,8 @@ function test_get_adjacent_to_quadratic_constraint() end -function test_get_adjacent_to_nonlinear_constraint() - m = make_degenerate_flow_model() +function test_get_adjacent_to_nonlinear_constraint(model_function=make_degenerate_flow_model) + m = model_function() igraph = MathProgIncidence.IncidenceGraphInterface(m) con = m[:bulk_dens_eqn] adjacent = MathProgIncidence.get_adjacent(igraph, con) @@ -135,8 +135,8 @@ function test_get_adjacent_to_nonlinear_constraint() end -function test_get_adjacent_to_variable() - m = make_degenerate_flow_model() +function test_get_adjacent_to_variable(model_function=make_degenerate_flow_model) + m = model_function() igraph = MathProgIncidence.IncidenceGraphInterface(m) var = m[:x][2] adjacent = MathProgIncidence.get_adjacent(igraph, var) @@ -151,8 +151,8 @@ function test_get_adjacent_to_variable() end -function test_maximum_matching() - m = make_degenerate_flow_model() +function test_maximum_matching(model_function=make_degenerate_flow_model) + m = model_function() igraph = MathProgIncidence.IncidenceGraphInterface(m) matching = MathProgIncidence.maximum_matching(igraph) @test length(matching) == 7 @@ -190,8 +190,8 @@ function test_maximum_matching() end -function test_dulmage_mendelsohn() - m = make_degenerate_flow_model() +function test_dulmage_mendelsohn(model_function=make_degenerate_flow_model) + m = model_function() igraph = MathProgIncidence.IncidenceGraphInterface(m) con_dmp, var_dmp = MathProgIncidence.dulmage_mendelsohn(igraph) con_undercon = con_dmp.underconstrained @@ -290,8 +290,8 @@ function test_dulmage_mendelsohn_from_constraints_and_variables() return end -function test_one_connected_component_igraph() - m = make_degenerate_flow_model() +function test_one_connected_component_igraph(model_function=make_degenerate_flow_model) + m = model_function() igraph = MathProgIncidence.IncidenceGraphInterface(m) con_comps, var_comps = MathProgIncidence.connected_components(igraph) @test length(var_comps) == 1 @@ -322,8 +322,8 @@ function test_multiple_connected_components_igraph() return end -function test_one_connected_component_cons_vars() - m = make_degenerate_flow_model() +function test_one_connected_component_cons_vars(model_function=make_degenerate_flow_model) + m = model_function() igraph = MathProgIncidence.IncidenceGraphInterface(m) con_dmp, var_dmp = MathProgIncidence.dulmage_mendelsohn(igraph) uc_var = [var_dmp.unmatched..., var_dmp.underconstrained...] @@ -351,8 +351,8 @@ function test_one_connected_component_cons_vars() return end -function test_construct_interface_active_inequalities() - m = make_simple_model() +function test_construct_interface_active_inequalities(model_function=make_simple_model) + m = model_function() JuMP.set_optimizer(m, Ipopt.Optimizer) JuMP.optimize!(m) @@ -375,8 +375,8 @@ function test_construct_interface_active_inequalities() return end -function test_active_inequalities_no_solution() - m = make_simple_model() +function test_active_inequalities_no_solution(model_function=make_simple_model) + m = model_function() @test_throws(JuMP.OptimizeNotCalled, igraph = MathProgIncidence.IncidenceGraphInterface( m, include_active_inequalities = true, tolerance = 1e-6 ) @@ -384,8 +384,8 @@ function test_active_inequalities_no_solution() return end -function test_bad_arguments() - m = make_simple_model() +function test_bad_arguments(model_function=make_simple_model) + m = model_function() @test_throws(ArgumentError, igraph = MathProgIncidence.IncidenceGraphInterface( m, include_active_inequalities = true, include_inequality = true ) @@ -395,22 +395,40 @@ end @testset "interface" begin test_construct_interface() + test_construct_interface(make_degenerate_flow_model_with_ScalarNonlinearFunction) test_construct_interface_rectangular() + test_construct_interface_rectangular(make_degenerate_flow_model_with_ScalarNonlinearFunction) test_get_adjacent_to_linear_constraint() + test_get_adjacent_to_linear_constraint(make_degenerate_flow_model_with_ScalarNonlinearFunction) test_get_adjacent_to_quadratic_constraint() + test_get_adjacent_to_quadratic_constraint(make_degenerate_flow_model_with_ScalarNonlinearFunction) test_get_adjacent_to_nonlinear_constraint() + test_get_adjacent_to_nonlinear_constraint(make_degenerate_flow_model_with_ScalarNonlinearFunction) test_get_adjacent_to_variable() + test_get_adjacent_to_variable(make_degenerate_flow_model_with_ScalarNonlinearFunction) test_maximum_matching() + test_maximum_matching(make_degenerate_flow_model_with_ScalarNonlinearFunction) test_dulmage_mendelsohn() + test_dulmage_mendelsohn(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_overconstrained_due_to_fixed_variable() test_overconstrained_due_to_including_bound() test_interface_from_constraints_and_variables() test_matching_from_constraints_and_variables() test_dulmage_mendelsohn_from_constraints_and_variables() - test_one_connected_component_igraph() test_multiple_connected_components_igraph() + + test_one_connected_component_igraph() + test_one_connected_component_igraph(make_degenerate_flow_model_with_ScalarNonlinearFunction) test_one_connected_component_cons_vars() + test_one_connected_component_cons_vars(make_degenerate_flow_model_with_ScalarNonlinearFunction) + test_construct_interface_active_inequalities() + test_construct_interface_active_inequalities(make_simple_model_with_ScalarNonlinearFunction) + test_active_inequalities_no_solution() + test_active_inequalities_no_solution(make_simple_model_with_ScalarNonlinearFunction) + test_bad_arguments() + test_bad_arguments(make_simple_model_with_ScalarNonlinearFunction) end diff --git a/test/models.jl b/test/models.jl index 5d75045..20bf1b4 100644 --- a/test/models.jl +++ b/test/models.jl @@ -38,6 +38,25 @@ function make_degenerate_flow_model() return m end +function make_degenerate_flow_model_with_ScalarNonlinearFunction() + m = JuMP.Model() + comps = [1, 2, 3] + @JuMP.variable(m, x[comps], start=1/3.0) + @JuMP.variable(m, flow_comp[comps], start=10.0) + @JuMP.variable(m, flow, start=30.0) + @JuMP.variable(m, rho, start=1.0) + + # sum_component_eqn + @JuMP.constraint(m, sum_comp_eqn, sum(x) == 1) + # component_density_eqn + @JuMP.constraint(m, comp_dens_eqn, x*rho .== [1.0, 1.1, 1.2]) + # density_eqn + @JuMP.constraint(m, bulk_dens_eqn, 1/rho - sum(1/x[j] for j in comps) == 0) + # component_flow_eqn + @JuMP.constraint(m, comp_flow_eqn, x.*flow .== flow_comp) + return m +end + function make_simple_model() m = JuMP.Model() JuMP.@variable(m, x[1:3] >= 0, start = 1.0) @@ -48,3 +67,14 @@ function make_simple_model() JuMP.@constraint(m, range1, 0 <= x[2] + x[3] <= 10) return m end + +function make_simple_model_with_ScalarNonlinearFunction() + m = JuMP.Model() + JuMP.@variable(m, x[1:3] >= 0, start = 1.0) + JuMP.@objective(m, Min, x[1]^2 + 2*x[2]^2 + 3*x[3]^2) + JuMP.@constraint(m, eq1, x[1]^1.5 * x[2]^2 == x[3]) + JuMP.@constraint(m, ineq1, - x[3] <= - 1.0) + JuMP.@constraint(m, ineq2, x[1] + 2*x[2] >= 4) + JuMP.@constraint(m, range1, 0 <= x[2] + x[3] <= 10) + return m +end