From 868dc27e6941a2126fcd09807bc10c2f44aabc39 Mon Sep 17 00:00:00 2001 From: odow Date: Fri, 30 Apr 2021 16:52:06 +1200 Subject: [PATCH 01/12] Add NL file writer. This is a single commit that copies the NL file writer from AmplNLWriter.jl that I recently wrote. It's a little different from the other FileFormats models in that it only supports copy_to, and it doesn't support reading. --- docs/src/submodules/FileFormats/overview.md | 7 + src/FileFormats/FileFormats.jl | 8 +- src/FileFormats/NL/NL.jl | 796 ++++++++++++++++++++ src/FileFormats/NL/NLExpr.jl | 492 ++++++++++++ src/FileFormats/NL/opcode.jl | 72 ++ test/FileFormats/FileFormats.jl | 4 +- test/FileFormats/NL/NL.jl | 759 +++++++++++++++++++ 7 files changed, 2136 insertions(+), 2 deletions(-) create mode 100644 src/FileFormats/NL/NL.jl create mode 100644 src/FileFormats/NL/NLExpr.jl create mode 100644 src/FileFormats/NL/opcode.jl create mode 100644 test/FileFormats/NL/NL.jl diff --git a/docs/src/submodules/FileFormats/overview.md b/docs/src/submodules/FileFormats/overview.md index 1c9c062ba9..b7d77c4318 100644 --- a/docs/src/submodules/FileFormats/overview.md +++ b/docs/src/submodules/FileFormats/overview.md @@ -45,6 +45,13 @@ julia> model = MOI.FileFormats.Model(format = MOI.FileFormats.FORMAT_MPS) A Mathematical Programming System (MPS) model ``` +**The NL file format** + +```jldoctest +julia> model = MOI.FileFormats.Model(format = MOI.FileFormats.FORMAT_NL) +An AMPL (.nl) model +``` + **The SDPA file format** ```jldoctest diff --git a/src/FileFormats/FileFormats.jl b/src/FileFormats/FileFormats.jl index 295bccc3f4..7758d81a13 100644 --- a/src/FileFormats/FileFormats.jl +++ b/src/FileFormats/FileFormats.jl @@ -12,6 +12,7 @@ include("CBF/CBF.jl") include("LP/LP.jl") include("MOF/MOF.jl") include("MPS/MPS.jl") +include("NL/NL.jl") include("SDPA/SDPA.jl") """ @@ -24,6 +25,7 @@ List of accepted export formats. - `FORMAT_LP`: the LP file format - `FORMAT_MOF`: the MathOptFormat file format - `FORMAT_MPS`: the MPS file format +- `FORMAT_NL`: the AMPL .nl file format - `FORMAT_SDPA`: the SemiDefinite Programming Algorithm format """ @enum( @@ -33,6 +35,7 @@ List of accepted export formats. FORMAT_LP, FORMAT_MOF, FORMAT_MPS, + FORMAT_NL, FORMAT_SDPA, ) @@ -64,6 +67,8 @@ function Model(; return MOF.Model(; kwargs...) elseif format == FORMAT_MPS return MPS.Model(; kwargs...) + elseif format == FORMAT_NL + return NL.Model(; kwargs...) elseif format == FORMAT_SDPA return SDPA.Model(; kwargs...) else @@ -76,6 +81,7 @@ function Model(; (".lp", LP.Model), (".mof.json", MOF.Model), (".mps", MPS.Model), + (".nl", NL.Model), (".dat-s", SDPA.Model), (".sdpa", SDPA.Model), ] @@ -88,7 +94,7 @@ function Model(; end const MATH_OPT_FORMATS = - Union{CBF.Model,LP.Model,MOF.Model,MPS.Model,SDPA.Model} + Union{CBF.Model,LP.Model,MOF.Model,MPS.Model,NL.Model,SDPA.Model} function MOI.write_to_file(model::MATH_OPT_FORMATS, filename::String) compressed_open(filename, "w", AutomaticCompression()) do io diff --git a/src/FileFormats/NL/NL.jl b/src/FileFormats/NL/NL.jl new file mode 100644 index 0000000000..0c77177adf --- /dev/null +++ b/src/FileFormats/NL/NL.jl @@ -0,0 +1,796 @@ +module NL + +import MathOptInterface +const MOI = MathOptInterface + +include("NLExpr.jl") + +### ============================================================================ +### Nonlinear constraints +### ============================================================================ + +struct _NLConstraint + lower::Float64 + upper::Float64 + opcode::Int + expr::_NLExpr +end + +""" + _NLConstraint(expr::Expr, bound::MOI.NLPBoundsPair) + +Convert a constraint in the form of a `expr` into a `_NLConstraint` object. + +See `MOI.constraint_expr` for details on the format. + +As a validation step, the right-hand side of each constraint must be a constant +term that is given by the `bound`. (If the constraint is an interval constraint, +both the left-hand and right-hand sides must be constants.) + +The six NL constraint types are: + + l <= g(x) <= u : 0 + g(x) >= l : 1 + g(x) <= u : 2 + g(x) : 3 # We don't support this + g(x) == c : 4 + x ⟂ g(x) : 5 # TODO(odow): Complementarity constraints +""" +function _NLConstraint(expr::Expr, bound::MOI.NLPBoundsPair) + if expr.head == :comparison + @assert length(expr.args) == 5 + if !(expr.args[1] ≈ bound.lower && bound.upper ≈ expr.args[5]) + _warn_invalid_bound(expr, bound) + end + return _NLConstraint( + expr.args[1], + expr.args[5], + 0, + _NLExpr(expr.args[3]), + ) + else + @assert expr.head == :call + @assert length(expr.args) == 3 + if expr.args[1] == :(<=) + if !(-Inf ≈ bound.lower && bound.upper ≈ expr.args[3]) + _warn_invalid_bound(expr, bound) + end + return _NLConstraint(-Inf, expr.args[3], 1, _NLExpr(expr.args[2])) + elseif expr.args[1] == :(>=) + if !(expr.args[3] ≈ bound.lower && bound.upper ≈ Inf) + _warn_invalid_bound(expr, bound) + end + return _NLConstraint(expr.args[3], Inf, 2, _NLExpr(expr.args[2])) + else + @assert expr.args[1] == :(==) + if !(expr.args[3] ≈ bound.lower ≈ bound.upper) + _warn_invalid_bound(expr, bound) + end + return _NLConstraint( + expr.args[3], + expr.args[3], + 4, + _NLExpr(expr.args[2]), + ) + end + end +end + +function _warn_invalid_bound(expr::Expr, bound::MOI.NLPBoundsPair) + return @warn( + "Invalid bounds detected in nonlinear constraint. Expected " * + "`$(bound.lower) <= g(x) <= $(bound.upper)`, but got the constraint " * + "$(expr)", + ) +end + +### ============================================================================ +### Nonlinear models +### ============================================================================ + +@enum(_VariableType, _BINARY, _INTEGER, _CONTINUOUS) + +mutable struct _VariableInfo + # Variable lower bound. + lower::Float64 + # Variable upper bound. + upper::Float64 + # Whether variable is binary or integer. + type::_VariableType + # Primal start of the variable. + start::Union{Float64,Nothing} + # Number of constraints that the variable appears in. + jacobian_count::Int + # If the variable appears in the objective. + in_nonlinear_objective::Bool + # If the objetive appears in a nonlinear constraint. + in_nonlinear_constraint::Bool + # The 0-indexed column of the variable. Computed right at the end. + order::Int + function _VariableInfo() + return new(-Inf, Inf, _CONTINUOUS, nothing, 0, false, false, 0) + end +end + +mutable struct Model <: MOI.ModelLike + # Store MOI.Name(). + name::String + # The objective expression. + f::_NLExpr + sense::MOI.OptimizationSense + # Number of nonlinear constraints in NLPBlock + nlpblock_dim::Int + # A vector of nonlinear constraints + g::Vector{_NLConstraint} + # A vector of linear constraints + h::Vector{_NLConstraint} + # A dictionary of info for the variables. + x::Dict{MOI.VariableIndex,_VariableInfo} + # A struct to help sort the mess that is variable ordering in NL files. + types::Vector{Vector{MOI.VariableIndex}} + # A vector of the final ordering of the variables. + order::Vector{MOI.VariableIndex} + + """ + Model() + + Create a new Optimizer object. + """ + function Model(; kwargs...) + return new( + "", + _NLExpr(false, _NLTerm[], Dict{MOI.VariableIndex,Float64}(), 0.0), + MOI.FEASIBILITY_SENSE, + 0, + _NLConstraint[], + _NLConstraint[], + Dict{MOI.VariableIndex,_VariableInfo}(), + [MOI.VariableIndex[] for _ in 1:9], + MOI.VariableIndex[], + ) + end +end + +Base.show(io::IO, ::Model) = print(io, "An AMPL (.nl) model") + +MOI.get(model::Model, ::MOI.SolverName) = "AmplNLWriter" + +MOI.supports(::Model, ::MOI.NLPBlock) = true + +MOI.supports(::Model, ::MOI.Name) = true +MOI.get(model::Model, ::MOI.Name) = model.name +MOI.set(model::Model, ::MOI.Name, name::String) = (model.name = name) + +function MOI.empty!(model::Model) + model.results = _NLResults( + "Optimize not called.", + MOI.OPTIMIZE_NOT_CALLED, + MOI.NO_SOLUTION, + NaN, + Dict{MOI.VariableIndex,Float64}(), + Float64[], + Dict{MOI.VariableIndex,Float64}(), + Dict{MOI.VariableIndex,Float64}(), + ) + model.f = _NLExpr(false, _NLTerm[], Dict{MOI.VariableIndex,Float64}(), 0.0) + empty!(model.g) + model.nlpblock_dim = 0 + empty!(model.h) + empty!(model.x) + for i in 1:9 + empty!(model.types[i]) + end + empty!(model.order) + return +end + +function MOI.is_empty(model::Model) + return isempty(model.g) && isempty(model.h) && isempty(model.x) +end + +const _SCALAR_FUNCTIONS = Union{ + MOI.SingleVariable, + MOI.ScalarAffineFunction{Float64}, + MOI.ScalarQuadraticFunction{Float64}, +} + +const _SCALAR_SETS = Union{ + MOI.LessThan{Float64}, + MOI.GreaterThan{Float64}, + MOI.EqualTo{Float64}, + MOI.Interval{Float64}, +} + +function MOI.supports_constraint( + ::Model, + ::Type{<:_SCALAR_FUNCTIONS}, + ::Type{<:_SCALAR_SETS}, +) + return true +end + +function MOI.supports_constraint( + ::Model, + ::Type{MOI.SingleVariable}, + ::Type{<:Union{MOI.ZeroOne,MOI.Integer}}, +) + return true +end + +MOI.supports(::Model, ::MOI.ObjectiveSense) = true +MOI.supports(::Model, ::MOI.ObjectiveFunction{<:_SCALAR_FUNCTIONS}) = true + +MOI.supports(::Model, ::MOI.RawParameter) = true +function MOI.get(model::Model, attr::MOI.RawParameter) + return get(model.options, attr.name, nothing) +end +function MOI.set(model::Model, attr::MOI.RawParameter, value) + model.options[attr.name] = value + return +end + +# ============================================================================== + +function MOI.supports( + ::Model, + ::MOI.VariablePrimalStart, + ::Type{MOI.VariableIndex}, +) + return true +end + +function MOI.set(model::Model, ::MOI.VariablePrimalStart, x, v::Real) + model.x[x].start = Float64(v) + return +end + +function MOI.set(model::Model, ::MOI.VariablePrimalStart, x, ::Nothing) + model.x[x].start = nothing + return +end + +function MOI.get(model::Model, ::MOI.VariablePrimalStart, x::MOI.VariableIndex) + return model.x[x].start +end + +# ============================================================================== + +struct _LinearNLPEvaluator <: MOI.AbstractNLPEvaluator end +MOI.features_available(::_LinearNLPEvaluator) = [:ExprGraph] +MOI.initialize(::_LinearNLPEvaluator, ::Vector{Symbol}) = nothing + +MOI.Utilities.supports_default_copy_to(::Model, ::Bool) = false + +function MOI.copy_to( + dest::Model, + model::MOI.ModelLike; + copy_names::Bool = false, +) + mapping = MOI.Utilities.IndexMap() + # Initialize the NLP block. + nlp_block = try + MOI.get(model, MOI.NLPBlock()) + catch + MOI.NLPBlockData(MOI.NLPBoundsPair[], _LinearNLPEvaluator(), false) + end + if !(:ExprGraph in MOI.features_available(nlp_block.evaluator)) + error( + "Unable to use AmplNLWriter because the nonlinear evaluator " * + "does not supply expression graphs.", + ) + end + MOI.initialize(nlp_block.evaluator, [:ExprGraph]) + # Objective function. + if nlp_block.has_objective + dest.f = _NLExpr(MOI.objective_expr(nlp_block.evaluator)) + else + F = MOI.get(model, MOI.ObjectiveFunctionType()) + obj = MOI.get(model, MOI.ObjectiveFunction{F}()) + dest.f = _NLExpr(obj) + end + # Nonlinear constraints + for (i, bound) in enumerate(nlp_block.constraint_bounds) + push!( + dest.g, + _NLConstraint(MOI.constraint_expr(nlp_block.evaluator, i), bound), + ) + end + dest.nlpblock_dim = length(dest.g) + for x in MOI.get(model, MOI.ListOfVariableIndices()) + dest.x[x] = _VariableInfo() + start = MOI.get(model, MOI.VariablePrimalStart(), x) + MOI.set(dest, MOI.VariablePrimalStart(), x, start) + mapping[x] = x + end + dest.sense = MOI.get(model, MOI.ObjectiveSense()) + resize!(dest.order, length(dest.x)) + # Now deal with the normal MOI constraints. + for (F, S) in MOI.get(model, MOI.ListOfConstraints()) + _process_constraint(dest, model, F, S, mapping) + end + # Correct bounds of binary variables. Mainly because AMPL doesn't have the + # concept of binary nonlinear variables, but it does have binary linear + # variables! How annoying. + for (x, v) in dest.x + if v.type == _BINARY + v.lower = max(0.0, v.lower) + v.upper = min(1.0, v.upper) + end + end + # Jacobian counts. The zero terms for nonlinear constraints should have + # been added when the expression was constructed. + for g in dest.g, v in keys(g.expr.linear_terms) + dest.x[v].jacobian_count += 1 + end + for h in dest.h, v in keys(h.expr.linear_terms) + dest.x[v].jacobian_count += 1 + end + # Now comes the confusing part. + # + # AMPL, in all its wisdom, orders variables in a _very_ specific way. + # The only hint in "Writing NL files" is the line "Variables are ordered as + # described in Tables 3 and 4 of [5]". + # + # Reading these + # + # https://cfwebprod.sandia.gov/cfdocs/CompResearch/docs/nlwrite20051130.pdf + # https://ampl.com/REFS/hooking2.pdf + # + # leads us to the following order + # + # 1) Continuous variables that appear in a + # nonlinear objective AND a nonlinear constraint + # 2) Discrete variables that appear in a + # nonlinear objective AND a nonlinear constraint + # 3) Continuous variables that appear in a + # nonlinear constraint, but NOT a nonlinear objective + # 4) Discrete variables that appear in a + # nonlinear constraint, but NOT a nonlinear objective + # 5) Continuous variables that appear in a + # nonlinear objective, but NOT a nonlinear constraint + # 6) Discrete variables that appear in a + # nonlinear objective, but NOT a nonlinear constraint + # 7) Continuous variables that DO NOT appear in a + # nonlinear objective or a nonlinear constraint + # 8) Binary variables that DO NOT appear in a + # nonlinear objective or a nonlinear constraint + # 9) Integer variables that DO NOT appear in a + # nonlinear objective or a nonlinear constraint + # + # Yes, nonlinear variables are broken into continuous/discrete, but linear + # variables are partitioned into continuous, binary, and integer. (See also, + # the need to modify bounds for binary variables.) + if !dest.f.is_linear + for x in keys(dest.f.linear_terms) + dest.x[x].in_nonlinear_objective = true + end + for x in dest.f.nonlinear_terms + if x isa MOI.VariableIndex + dest.x[x].in_nonlinear_objective = true + end + end + end + for con in dest.g + for x in keys(con.expr.linear_terms) + dest.x[x].in_nonlinear_constraint = true + end + for x in con.expr.nonlinear_terms + if x isa MOI.VariableIndex + dest.x[x].in_nonlinear_constraint = true + end + end + end + types = dest.types + for (x, v) in dest.x + if v.in_nonlinear_constraint && v.in_nonlinear_objective + push!(v.type == _CONTINUOUS ? types[1] : types[2], x) + elseif v.in_nonlinear_constraint + push!(v.type == _CONTINUOUS ? types[3] : types[4], x) + elseif v.in_nonlinear_objective + push!(v.type == _CONTINUOUS ? types[5] : types[6], x) + elseif v.type == _CONTINUOUS + push!(types[7], x) + elseif v.type == _BINARY + push!(types[8], x) + else + @assert v.type == _INTEGER + push!(types[9], x) + end + end + # However! Don't let Tables 3 and 4 fool you, because the ordering actually + # depends on whether the number of nonlinear variables in the objective only + # is _strictly_ greater than the number of nonlinear variables in the + # constraints only. Quoting: + # + # For all versions, the first nlvc variables appear nonlinearly in at + # least one constraint. If nlvo > nlvc, the first nlvc variables may or + # may not appear nonlinearly in an objective, but the next nlvo – nlvc + # variables do appear nonlinearly in at least one objective. Otherwise + # all of the first nlvo variables appear nonlinearly in an objective. + # + # However, even this is slightly incorrect, because I think it should read + # "For all versions, the first nlvb variables appear nonlinearly." The "nlvo + # - nlvc" part is also clearly incorrect, and should probably read "nlvo - + # nlvb." + # + # It's a bit confusing, so here is the relevant code from Couenne: + # https://github.com/coin-or/Couenne/blob/683c5b305d78a009d59268a4bca01e0ad75ebf02/src/readnl/readnl.cpp#L76-L87 + # + # They interpret this paragraph to mean the switch on nlvo > nlvc determines + # whether the next block of variables are the ones that appear in the + # objective only, or the constraints only. + # + # That makes sense as a design choice, because you can read them in two + # contiguous blocks. + # + # Essentially, what all this means is if !(nlvo > nlvc), then swap 3-4 for + # 5-6 in the variable order. + nlvc = length(types[3]) + length(types[4]) + nlvo = length(types[5]) + length(types[6]) + order_i = if nlvo > nlvc + [1, 2, 3, 4, 5, 6, 7, 8, 9] + else + [1, 2, 5, 6, 3, 4, 7, 8, 9] + end + # Now we can order the variables. + n = 0 + for i in order_i + # Since variables come from a dictionary, there may be differences in + # the order depending on platform and Julia version. Sort by creation + # time for consistency. + for x in sort!(types[i]; by = y -> y.value) + dest.x[x].order = n + dest.order[n+1] = x + n += 1 + end + end + return mapping +end + +_set_to_bounds(set::MOI.Interval) = (0, set.lower, set.upper) +_set_to_bounds(set::MOI.LessThan) = (1, -Inf, set.upper) +_set_to_bounds(set::MOI.GreaterThan) = (2, set.lower, Inf) +_set_to_bounds(set::MOI.EqualTo) = (4, set.value, set.value) + +function _process_constraint(dest::Model, model, F, S, mapping) + for ci in MOI.get(model, MOI.ListOfConstraintIndices{F,S}()) + f = MOI.get(model, MOI.ConstraintFunction(), ci) + s = MOI.get(model, MOI.ConstraintSet(), ci) + op, l, u = _set_to_bounds(s) + con = _NLConstraint(l, u, op, _NLExpr(f)) + if isempty(con.expr.linear_terms) && isempty(con.expr.nonlinear_terms) + if l <= con.expr.constant <= u + continue + else + error( + "Malformed constraint. There are no variables and the " * + "bounds don't make sense.", + ) + end + elseif con.expr.is_linear + push!(dest.h, con) + mapping[ci] = MOI.ConstraintIndex{F,S}(length(dest.h)) + else + push!(dest.g, con) + mapping[ci] = MOI.ConstraintIndex{F,S}(length(dest.g)) + end + end + return +end + +function _process_constraint( + dest::Model, + model, + F::Type{MOI.SingleVariable}, + S, + mapping, +) + for ci in MOI.get(model, MOI.ListOfConstraintIndices{F,S}()) + mapping[ci] = ci + f = MOI.get(model, MOI.ConstraintFunction(), ci) + s = MOI.get(model, MOI.ConstraintSet(), ci) + _, l, u = _set_to_bounds(s) + if l > -Inf + dest.x[f.variable].lower = l + end + if u < Inf + dest.x[f.variable].upper = u + end + end + return +end + +function _process_constraint( + dest::Model, + model, + F::Type{MOI.SingleVariable}, + S::Union{Type{MOI.ZeroOne},Type{MOI.Integer}}, + mapping, +) + for ci in MOI.get(model, MOI.ListOfConstraintIndices{F,S}()) + mapping[ci] = ci + f = MOI.get(model, MOI.ConstraintFunction(), ci) + dest.x[f.variable].type = S == MOI.ZeroOne ? _BINARY : _INTEGER + end + return +end + +_str(x::Float64) = isinteger(x) ? string(round(Int, x)) : string(x) + +_write_term(io, x::Float64, ::Any) = println(io, "n", _str(x)) +_write_term(io, x::Int, ::Any) = println(io, "o", x) +function _write_term(io, x::MOI.VariableIndex, nlmodel) + return println(io, "v", nlmodel.x[x].order) +end + +_is_nary(x::Int) = x in _NARY_OPCODES +_is_nary(x) = false + +function _write_nlexpr(io::IO, expr::_NLExpr, nlmodel::Model) + if expr.is_linear || length(expr.nonlinear_terms) == 0 + # If the expression is linear, just write out the constant term. + _write_term(io, expr.constant, nlmodel) + return + end + # If the constant term is non-zero, we need to write it out. + skip_terms = 0 + if !iszero(expr.constant) + if expr.nonlinear_terms[1] == OPSUMLIST + # The nonlinear expression is a summation. We can write our constant + # first, but we also need to increment the number of arguments by + # one. In addition, since we're writing out the first two terms now, + # we must skip them below. + _write_term(io, OPSUMLIST, nlmodel) + println(io, expr.nonlinear_terms[2] + 1) + _write_term(io, expr.constant, nlmodel) + skip_terms = 2 + else + # The nonlinear expression is something other than a summation, so + # add a new + node to the expression. + _write_term(io, OPPLUS, nlmodel) + _write_term(io, expr.constant, nlmodel) + end + end + last_nary = false + for term in expr.nonlinear_terms + if skip_terms > 0 + skip_terms -= 1 + continue + end + if last_nary + println(io, term::Int) + last_nary = false + else + _write_term(io, term, nlmodel) + last_nary = _is_nary(term) + end + end + return +end + +function _write_linear_block(io::IO, expr::_NLExpr, nlmodel::Model) + elements = [(c, nlmodel.x[v].order) for (v, c) in expr.linear_terms] + for (c, x) in sort!(elements; by = i -> i[2]) + println(io, x, " ", _str(c)) + end + return +end + +function Base.write(io::IO, nlmodel::Model) + # ========================================================================== + # Header + # Line 1: Always the same + # Notes: + # * I think these are magic bytes used by AMPL internally for stuff. + println(io, "g3 1 1 0") + + # Line 2: vars, constraints, objectives, ranges, eqns, logical constraints + # Notes: + # * We assume there is always one objective, even if it is just `min 0`. + n_con, n_ranges, n_eqns = 0, 0, 0 + for cons in (nlmodel.g, nlmodel.h), c in cons + n_con += 1 + if c.opcode == 0 + n_ranges += 1 + elseif c.opcode == 4 + n_eqns += 1 + end + end + println(io, " $(length(nlmodel.x)) $(n_con) 1 $(n_ranges) $(n_eqns) 0") + + # Line 3: nonlinear constraints, objectives + # Notes: + # * We assume there is always one objective, even if it is just `min 0`. + n_nlcon = length(nlmodel.g) + println(io, " ", n_nlcon, " ", 1) + + # Line 4: network constraints: nonlinear, linear + # Notes: + # * We don't support network constraints. I don't know how they are + # represented. + println(io, " 0 0") + + # Line 5: nonlinear vars in constraints, objectives, both + # Notes: + # * This order is confusingly different to the standard "b, c, o" order. + nlvb = length(nlmodel.types[1]) + length(nlmodel.types[2]) + nlvc = nlvb + length(nlmodel.types[3]) + length(nlmodel.types[4]) + nlvo = nlvb + length(nlmodel.types[5]) + length(nlmodel.types[6]) + println(io, " ", nlvc, " ", nlvo, " ", nlvb) + + # Line 6: linear network variables; functions; arith, flags + # Notes: + # * I don't know what this line means. It is what it is. Apparently `flags` + # is set to 1 to get suffixes in .sol file. + println(io, " 0 0 0 1") + + # Line 7: discrete variables: binary, integer, nonlinear (b,c,o) + # Notes: + # * The order is + # - binary variables in linear only + # - integer variables in linear only + # - binary or integer variables in nonlinear objective and constraint + # - binary or integer variables in nonlinear constraint + # - binary or integer variables in nonlinear objective + nbv = length(nlmodel.types[8]) + niv = length(nlmodel.types[9]) + nl_both = length(nlmodel.types[2]) + nl_cons = length(nlmodel.types[4]) + nl_obj = length(nlmodel.types[6]) + println(io, " ", nbv, " ", niv, " ", nl_both, " ", nl_cons, " ", nl_obj) + + # Line 8: nonzeros in Jacobian, gradients + # Notes: + # * Make sure to include a 0 element for every variable that appears in an + # objective or constraint, even if the linear coefficient is 0. + nnz_jacobian = 0 + for g in nlmodel.g + nnz_jacobian += length(g.expr.linear_terms) + end + for h in nlmodel.h + nnz_jacobian += length(h.expr.linear_terms) + end + nnz_gradient = length(nlmodel.f.linear_terms) + println(io, " ", nnz_jacobian, " ", nnz_gradient) + + # Line 9: max name lengths: constraints, variables + # Notes: + # * We don't add names, so this is just 0, 0. + println(io, " 0 0") + + # Line 10: common exprs: b,c,o,c1,o1 + # Notes: + # * We don't add common subexpressions (i.e., V blocks). + # * I assume the notation means + # - b = in nonlinear objective and constraint + # - c = in nonlinear constraint + # - o = in nonlinear objective + # - c1 = in linear constraint + # - o1 = in linear objective + println(io, " 0 0 0 0 0") + # ========================================================================== + # Constraints + # Notes: + # * Nonlinear constraints first, then linear. + # * For linear constraints, write out the constant term here. + for (i, g) in enumerate(nlmodel.g) + println(io, "C", i - 1) + _write_nlexpr(io, g.expr, nlmodel) + end + for (i, h) in enumerate(nlmodel.h) + println(io, "C", i - 1 + n_nlcon) + _write_nlexpr(io, h.expr, nlmodel) + end + # ========================================================================== + # Objective + # Notes: + # * NL files support multiple objectives, but we're just going to write 1, + # so it's always `O0`. + # * For linear objectives, write out the constant term here. + println(io, "O0 ", nlmodel.sense == MOI.MAX_SENSE ? "1" : "0") + _write_nlexpr(io, nlmodel.f, nlmodel) + # ========================================================================== + # VariablePrimalStart + # Notes: + # * Make sure to write out the variables in order. + println(io, "x", length(nlmodel.x)) + for (i, x) in enumerate(nlmodel.order) + start = MOI.get(nlmodel, MOI.VariablePrimalStart(), x) + println(io, i - 1, " ", start === nothing ? 0 : _str(start)) + end + # ========================================================================== + # Constraint bounds + # Notes: + # * Nonlinear constraints go first, then linear. + # * The constant term for linear constraints gets written out in the + # "C" block. + if n_con > 0 + println(io, "r") + # Nonlinear constraints + for g in nlmodel.g + print(io, g.opcode) + if g.opcode == 0 + println(io, " ", _str(g.lower), " ", _str(g.upper)) + elseif g.opcode == 1 + println(io, " ", _str(g.upper)) + elseif g.opcode == 2 + println(io, " ", _str(g.lower)) + else + @assert g.opcode == 4 + println(io, " ", _str(g.lower)) + end + end + # Linear constraints + for h in nlmodel.h + print(io, h.opcode) + if h.opcode == 0 + println(io, " ", _str(h.lower), " ", _str(h.upper)) + elseif h.opcode == 1 + println(io, " ", _str(h.upper)) + elseif h.opcode == 2 + println(io, " ", _str(h.lower)) + else + @assert h.opcode == 4 + println(io, " ", _str(h.lower)) + end + end + end + # ========================================================================== + # Variable bounds + # Notes: + # * Not much to note, other than to make sure you iterate the variables in + # the correct order. + println(io, "b") + for x in nlmodel.order + v = nlmodel.x[x] + if v.lower == v.upper + println(io, "4 ", _str(v.lower)) + elseif -Inf < v.lower && v.upper < Inf + println(io, "0 ", _str(v.lower), " ", _str(v.upper)) + elseif -Inf == v.lower && v.upper < Inf + println(io, "1 ", _str(v.upper)) + elseif -Inf < v.lower && v.upper == Inf + println(io, "2 ", _str(v.lower)) + else + println(io, "3") + end + end + # ========================================================================== + # Jacobian block + # Notes: + # * If a variable appears in a constraint, it needs to have a corresponding + # entry in the Jacobian block, even if the linear coefficient is zero. + # AMPL uses this to determine the Jacobian sparsity. + # * As before, nonlinear constraints go first, then linear. + # * Don't write out the `k` entry for the last variable, because it can be + # inferred from the total number of elements in the Jacobian as given in + # the header. + if n_con > 0 + println(io, "k", length(nlmodel.x) - 1) + total = 0 + for i in 1:length(nlmodel.order)-1 + total += nlmodel.x[nlmodel.order[i]].jacobian_count + println(io, total) + end + for (i, g) in enumerate(nlmodel.g) + println(io, "J", i - 1, " ", length(g.expr.linear_terms)) + _write_linear_block(io, g.expr, nlmodel) + end + for (i, h) in enumerate(nlmodel.h) + println(io, "J", i - 1 + n_nlcon, " ", length(h.expr.linear_terms)) + _write_linear_block(io, h.expr, nlmodel) + end + end + # ========================================================================== + # Gradient block + # Notes: + # * You only need to write this out if there are linear terms in the + # objective. + if nnz_gradient > 0 + println(io, "G0 ", nnz_gradient) + _write_linear_block(io, nlmodel.f, nlmodel) + end + return nlmodel +end + +end diff --git a/src/FileFormats/NL/NLExpr.jl b/src/FileFormats/NL/NLExpr.jl new file mode 100644 index 0000000000..5dc39e3371 --- /dev/null +++ b/src/FileFormats/NL/NLExpr.jl @@ -0,0 +1,492 @@ +### ============================================================================ +### Opcodes, and their AMPL <-> Julia conversions. +### ============================================================================ + +include("opcode.jl") + +""" + _JULIA_TO_AMPL + +This dictionary is manualy curated, based on the list of opcodes in `opcode.jl`. + +The goal is to map Julia functions to their AMPL opcode equivalent. + +Sometimes, there is ambiguity, such as the `:+`, which Julia uses for unary, +binary, and n-ary addition, while AMPL doesn't support unary addition, uses +OPPLUS for binary, and OPSUMLIST for n-ary. In these cases, introduce a +different symbol to disambiguate them in the context of this dictionary, and add +logic to `_process_expr!` to rewrite the Julia expression. + +Commented out lines are opcodes supported by AMPL that don't have a clear Julia +equivalent. If you can think of one, feel free to add it. But then go and make +similar changes to `_AMPL_TO_JULIA` and `_NARY_OPCODES`. +""" +const _JULIA_TO_AMPL = Dict{Symbol,Int}( + :+ => OPPLUS, # binary-plus + :- => OPMINUS, + :* => OPMULT, + :/ => OPDIV, + :rem => OPREM, + :^ => OPPOW, + # OPLESS = 6 + :min => MINLIST, # n-ary + :max => MAXLIST, # n-ary + # FLOOR = 13 + # CEIL = 14 + :abs => ABS, + :neg => OPUMINUS, + :|| => OPOR, + :&& => OPAND, + :(<) => LT, + :(<=) => LE, + :(==) => EQ, + :(>=) => GE, + :(>) => GT, + :(!=) => NE, + :(!) => OPNOT, + :ifelse => OPIFnl, + :tanh => OP_tanh, + :tan => OP_tan, + :sqrt => OP_sqrt, + :sinh => OP_sinh, + :sin => OP_sin, + :log10 => OP_log10, + :log => OP_log, + :exp => OP_exp, + :cosh => OP_cosh, + :cos => OP_cos, + :atanh => OP_atanh, + # OP_atan2 = 48, + :atan => OP_atan, + :asinh => OP_asinh, + :asin => OP_asin, + :acosh => OP_acosh, + :acos => OP_acos, + :sum => OPSUMLIST, # n-ary plus + # OPintDIV = 55 + # OPprecision = 56 + # OPround = 57 + # OPtrunc = 58 + # OPCOUNT = 59 + # OPNUMBEROF = 60 + # OPNUMBEROFs = 61 + # OPATLEAST = 62 + # OPATMOST = 63 + # OPPLTERM = 64 + # OPIFSYM = 65 + # OPEXACTLY = 66 + # OPNOTATLEAST = 67 + # OPNOTATMOST = 68 + # OPNOTEXACTLY = 69 + # ANDLIST = 70 + # ORLIST = 71 + # OPIMPELSE = 72 + # OP_IFF = 73 + # OPALLDIFF = 74 + # OPSOMESAME = 75 + # OP1POW = 76 + # OP2POW = 77 + # OPCPOW = 78 + # OPFUNCALL = 79 + # OPNUM = 80 + # OPHOL = 81 + # OPVARVAL = 82 + # N_OPS = 83 +) + +""" + _AMPL_TO_JULIA + +This dictionary is manualy curated, based on the list of supported opcodes +`_JULIA_TO_AMPL`. + +The goal is to map AMPL opcodes to their Julia equivalents. In addition, we need +to know the arity of each opcode. + +If the opcode is an n-ary function, use `-1`. +""" +const _AMPL_TO_JULIA = Dict{Int,Tuple{Int,Function}}( + OPPLUS => (2, +), + OPMINUS => (2, -), + OPMULT => (2, *), + OPDIV => (2, /), + OPREM => (2, rem), + OPPOW => (2, ^), + # OPLESS = 6 + MINLIST => (-1, minimum), + MAXLIST => (-1, maximum), + # FLOOR = 13 + # CEIL = 14 + ABS => (1, abs), + OPUMINUS => (1, -), + OPOR => (2, |), + OPAND => (2, &), + LT => (2, <), + LE => (2, <=), + EQ => (2, ==), + GE => (2, >=), + GT => (2, >), + NE => (2, !=), + OPNOT => (1, !), + OPIFnl => (3, ifelse), + OP_tanh => (1, tanh), + OP_tan => (1, tan), + OP_sqrt => (1, sqrt), + OP_sinh => (1, sinh), + OP_sin => (1, sin), + OP_log10 => (1, log10), + OP_log => (1, log), + OP_exp => (1, exp), + OP_cosh => (1, cosh), + OP_cos => (1, cos), + OP_atanh => (1, atanh), + # OP_atan2 = 48, + OP_atan => (1, atan), + OP_asinh => (1, asinh), + OP_asin => (1, asin), + OP_acosh => (1, acosh), + OP_acos => (1, acos), + OPSUMLIST => (-1, sum), + # OPintDIV = 55 + # OPprecision = 56 + # OPround = 57 + # OPtrunc = 58 + # OPCOUNT = 59 + # OPNUMBEROF = 60 + # OPNUMBEROFs = 61 + # OPATLEAST = 62 + # OPATMOST = 63 + # OPPLTERM = 64 + # OPIFSYM = 65 + # OPEXACTLY = 66 + # OPNOTATLEAST = 67 + # OPNOTATMOST = 68 + # OPNOTEXACTLY = 69 + # ANDLIST = 70 + # ORLIST = 71 + # OPIMPELSE = 72 + # OP_IFF = 73 + # OPALLDIFF = 74 + # OPSOMESAME = 75 + # OP1POW = 76 + # OP2POW = 77 + # OPCPOW = 78 + # OPFUNCALL = 79 + # OPNUM = 80 + # OPHOL = 81 + # OPVARVAL = 82 + # N_OPS = 83 +) + +""" + _NARY_OPCODES + +A manually curated list of n-ary opcodes, taken from Table 8 of "Writing .nl +files." + +Not all of these are implemented. See `_REV_OPCODES` for more detail. +""" +const _NARY_OPCODES = Set([ + MINLIST, + MAXLIST, + OPSUMLIST, + OPCOUNT, + OPNUMBEROF, + OPNUMBEROFs, + ANDLIST, + ORLIST, + OPALLDIFF, +]) + +""" + _UNARY_SPECIAL_CASES + +This dictionary defines a set of unary functions that are special-cased. They +don't exist in the NL file format, but they may be called from Julia, and +they can easily be converted into NL-compatible expressions. + +If you have a new unary-function that you want to support, add it here. +""" +const _UNARY_SPECIAL_CASES = Dict{Symbol,Function}( + :cbrt => (x) -> :($x^(1 / 3)), + :abs2 => (x) -> :($x^2), + :inv => (x) -> :(1 / $x), + :log2 => (x) -> :(log($x) / log(2)), + :log1p => (x) -> :(log(1 + $x)), + :exp2 => (x) -> :(2^$x), + :expm1 => (x) -> :(exp($x) - 1), + :sec => (x) -> :(1 / cos($x)), + :csc => (x) -> :(1 / sin($x)), + :cot => (x) -> :(1 / tan($x)), + :asec => (x) -> :(acos(1 / $x)), + :acsc => (x) -> :(asin(1 / $x)), + :acot => (x) -> :(pi / 2 - atan($x)), + :sind => (x) -> :(sin(pi / 180 * $x)), + :cosd => (x) -> :(cos(pi / 180 * $x)), + :tand => (x) -> :(tan(pi / 180 * $x)), + :secd => (x) -> :(1 / cos(pi / 180 * $x)), + :cscd => (x) -> :(1 / sin(pi / 180 * $x)), + :cotd => (x) -> :(1 / tan(pi / 180 * $x)), + :asind => (x) -> :(asin($x) * 180 / pi), + :acosd => (x) -> :(acos($x) * 180 / pi), + :atand => (x) -> :(atan($x) * 180 / pi), + :asecd => (x) -> :(acos(1 / $x) * 180 / pi), + :acscd => (x) -> :(asin(1 / $x) * 180 / pi), + :acotd => (x) -> :((pi / 2 - atan($x)) * 180 / pi), + :sech => (x) -> :(1 / cosh($x)), + :csch => (x) -> :(1 / sinh($x)), + :coth => (x) -> :(1 / tanh($x)), + :asech => (x) -> :(acosh(1 / $x)), + :acsch => (x) -> :(asinh(1 / $x)), + :acoth => (x) -> :(atanh(1 / $x)), +) + +""" + _BINARY_SPECIAL_CASES + +This dictionary defines a set of binary functions that are special-cased. They +don't exist in the NL file format, but they may be called from Julia, and +they can easily be converted into NL-compatible expressions. + +If you have a new binary-function that you want to support, add it here. +""" +const _BINARY_SPECIAL_CASES = Dict{Symbol,Function}(:\ => (x, y) -> :($y / $x)) + +### ============================================================================ +### Nonlinear expressions +### ============================================================================ + +# TODO(odow): This type isn't great. We should experiment with something that is +# type-stable, like +# +# @enum(_NLType, _INTEGER, _DOUBLE, _VARIABLE) +# struct _NLTerm +# type::_NLType +# data::Int64 +# end +# _NLTerm(x::Int) = _NLTerm(_INTEGER, x) +# _NLTerm(x::Float64) = _NLTerm(_DOUBLE, reinterpret(Int64, x)) +# _NLTerm(x::MOI.VariableIndex) = _NLTerm(_VARIABLE, x.value) +# function _value(x::_NLTerm) +# if x.type == _INTEGER +# return x.data +# elseif x.type == _DOUBLE +# return reinterpret(Float64, x.data) +# else +# @assert x.type == _VARIABLE +# return MOI.VariableIndex(x.data) +# end +# end + +const _NLTerm = Union{Int,Float64,MOI.VariableIndex} + +struct _NLExpr + is_linear::Bool + nonlinear_terms::Vector{_NLTerm} + linear_terms::Dict{MOI.VariableIndex,Float64} + constant::Float64 +end + +function Base.:(==)(x::_NLExpr, y::_NLExpr) + return x.is_linear == y.is_linear && + x.nonlinear_terms == y.nonlinear_terms && + x.linear_terms == y.linear_terms && + x.constant == y.constant +end + +_NLExpr(x::MOI.VariableIndex) = _NLExpr(true, _NLTerm[], Dict(x => 1.0), 0.0) + +_NLExpr(x::MOI.SingleVariable) = _NLExpr(x.variable) + +function _add_or_set(dict, key, value) + if haskey(dict, key) + dict[key] += value + else + dict[key] = value + end + return +end + +function _NLExpr(x::MOI.ScalarAffineFunction) + linear = Dict{MOI.VariableIndex,Float64}() + for (i, term) in enumerate(x.terms) + _add_or_set(linear, term.variable_index, term.coefficient) + end + return _NLExpr(true, _NLTerm[], linear, x.constant) +end + +function _NLExpr(x::MOI.ScalarQuadraticFunction) + linear = Dict{MOI.VariableIndex,Float64}() + for (i, term) in enumerate(x.affine_terms) + _add_or_set(linear, term.variable_index, term.coefficient) + end + terms = _NLTerm[] + N = length(x.quadratic_terms) + if N == 0 || N == 1 + # If there are 0 or 1 terms, no need for an addition node. + elseif N == 2 + # If there are two terms, use binary addition. + push!(terms, OPPLUS) + elseif N > 2 + # If there are more, use n-ary addition. + push!(terms, OPSUMLIST) + push!(terms, N) + end + for term in x.quadratic_terms + coefficient = term.coefficient + # MOI defines quadratic as 1/2 x' Q x :( + if term.variable_index_1 == term.variable_index_2 + coefficient *= 0.5 + end + # Optimization: no need for the OPMULT if the coefficient is 1. + if !isone(coefficient) + push!(terms, OPMULT) + push!(terms, coefficient) + end + push!(terms, OPMULT) + push!(terms, term.variable_index_1) + push!(terms, term.variable_index_2) + # For the Jacobian sparsity patterns, we need to add a linear term, even + # if the variable only appears nonlinearly. + _add_or_set(linear, term.variable_index_1, 0.0) + _add_or_set(linear, term.variable_index_2, 0.0) + end + return _NLExpr(false, terms, linear, x.constant) +end + +function _NLExpr(expr::Expr) + nlexpr = _NLExpr(false, _NLTerm[], Dict{MOI.VariableIndex,Float64}(), 0.0) + _process_expr!(nlexpr, expr) + return nlexpr +end + +function _process_expr!(expr::_NLExpr, arg::Real) + return push!(expr.nonlinear_terms, Float64(arg)) +end + +# Assume the symbol is a numeric constant like `pi`. +_process_expr!(expr::_NLExpr, arg::Symbol) = _process_expr!(expr, eval(arg)) + +function _process_expr!(expr::_NLExpr, arg::MOI.VariableIndex) + _add_or_set(expr.linear_terms, arg, 0.0) + return push!(expr.nonlinear_terms, arg) +end + +# TODO(odow): these process_expr! functions use recursion. For large models, +# this may exceed the stack. At some point, we may have to rewrite this to not +# use recursion. +function _process_expr!(expr::_NLExpr, arg::Expr) + if arg.head == :call + if length(arg.args) == 2 + f = get(_UNARY_SPECIAL_CASES, arg.args[1], nothing) + if f !== nothing + return _process_expr!(expr, f(arg.args[2])) + end + elseif length(arg.args) == 3 + f = get(_BINARY_SPECIAL_CASES, arg.args[1], nothing) + if f !== nothing + return _process_expr!(expr, f(arg.args[2], arg.args[3])) + end + end + return _process_expr!(expr, arg.args) + elseif arg.head == :ref + return _process_expr!(expr, arg.args[2]) + elseif arg == :() + return # Some evalators return a null objective of `:()`. + else + error("Unsupported expression: $(arg)") + end + return +end + +function _process_expr!(expr::_NLExpr, args::Vector{Any}) + op = first(args) + N = length(args) - 1 + # Before processing the arguments, do some re-writing. + if op == :+ + if N == 1 # +x, so we can just drop the op and process the args. + return _process_expr!(expr, args[2]) + elseif N > 2 # nary-addition! + op = :sum + end + elseif op == :- && N == 1 + op = :neg + elseif op == :* && N > 2 # nary-multiplication. + # NL doesn't define an nary multiplication operator, so we need to + # rewrite our expression as a sequence of chained binary operators. + while N > 2 + # Combine last term with previous to form a binary * expression + arg = pop!(args) + args[end] = Expr(:call, :*, args[end], arg) + N = length(args) - 1 + end + end + # Now convert the Julia expression into an _NLExpr. + opcode = get(_JULIA_TO_AMPL, op, nothing) + if opcode === nothing + error("Unsupported operation $(op)") + end + push!(expr.nonlinear_terms, opcode) + if opcode in _NARY_OPCODES + push!(expr.nonlinear_terms, N) + end + for i in 1:N + _process_expr!(expr, args[i+1]) + end + return +end + +### ============================================================================ +### Evaluate nonlinear expressions +### ============================================================================ + +function _evaluate(expr::_NLExpr, x::Dict{MOI.VariableIndex,Float64}) + y = expr.constant + for (v, c) in expr.linear_terms + y += c * x[v] + end + if length(expr.nonlinear_terms) > 0 + ret, n = _evaluate(expr.nonlinear_terms[1], expr.nonlinear_terms, x, 1) + @assert n == length(expr.nonlinear_terms) + 1 + y += ret + end + return y +end + +function _evaluate( + head::MOI.VariableIndex, + ::Vector{_NLTerm}, + x::Dict{MOI.VariableIndex,Float64}, + head_i::Int, +)::Tuple{Float64,Int} + return x[head], head_i + 1 +end + +function _evaluate( + head::Float64, + ::Vector{_NLTerm}, + ::Dict{MOI.VariableIndex,Float64}, + head_i::Int, +)::Tuple{Float64,Int} + return head, head_i + 1 +end + +function _evaluate( + head::Int, + terms::Vector{_NLTerm}, + x::Dict{MOI.VariableIndex,Float64}, + head_i::Int, +)::Tuple{Float64,Int} + N, f = _AMPL_TO_JULIA[head] + is_nary = (N == -1) + head_i += 1 + if is_nary + N = terms[head_i]::Int + head_i += 1 + end + args = Vector{Float64}(undef, N) + for n in 1:N + args[n], head_i = _evaluate(terms[head_i], terms, x, head_i) + end + return is_nary ? f(args) : f(args...), head_i +end diff --git a/src/FileFormats/NL/opcode.jl b/src/FileFormats/NL/opcode.jl new file mode 100644 index 0000000000..2a066b14ab --- /dev/null +++ b/src/FileFormats/NL/opcode.jl @@ -0,0 +1,72 @@ +# Do not modify. This file is automatically created by the script in `gen.jl`. + +const OPPLUS = 0 +const OPMINUS = 1 +const OPMULT = 2 +const OPDIV = 3 +const OPREM = 4 +const OPPOW = 5 +const OPLESS = 6 +const MINLIST = 11 +const MAXLIST = 12 +const FLOOR = 13 +const CEIL = 14 +const ABS = 15 +const OPUMINUS = 16 +const OPOR = 20 +const OPAND = 21 +const LT = 22 +const LE = 23 +const EQ = 24 +const GE = 28 +const GT = 29 +const NE = 30 +const OPNOT = 34 +const OPIFnl = 35 +const OP_tanh = 37 +const OP_tan = 38 +const OP_sqrt = 39 +const OP_sinh = 40 +const OP_sin = 41 +const OP_log10 = 42 +const OP_log = 43 +const OP_exp = 44 +const OP_cosh = 45 +const OP_cos = 46 +const OP_atanh = 47 +const OP_atan2 = 48 +const OP_atan = 49 +const OP_asinh = 50 +const OP_asin = 51 +const OP_acosh = 52 +const OP_acos = 53 +const OPSUMLIST = 54 +const OPintDIV = 55 +const OPprecision = 56 +const OPround = 57 +const OPtrunc = 58 +const OPCOUNT = 59 +const OPNUMBEROF = 60 +const OPNUMBEROFs = 61 +const OPATLEAST = 62 +const OPATMOST = 63 +const OPPLTERM = 64 +const OPIFSYM = 65 +const OPEXACTLY = 66 +const OPNOTATLEAST = 67 +const OPNOTATMOST = 68 +const OPNOTEXACTLY = 69 +const ANDLIST = 70 +const ORLIST = 71 +const OPIMPELSE = 72 +const OP_IFF = 73 +const OPALLDIFF = 74 +const OPSOMESAME = 75 +const OP1POW = 76 +const OP2POW = 77 +const OPCPOW = 78 +const OPFUNCALL = 79 +const OPNUM = 80 +const OPHOL = 81 +const OPVARVAL = 82 +const N_OPS = 83 diff --git a/test/FileFormats/FileFormats.jl b/test/FileFormats/FileFormats.jl index 08118db808..7e4d02afed 100644 --- a/test/FileFormats/FileFormats.jl +++ b/test/FileFormats/FileFormats.jl @@ -4,7 +4,7 @@ using Test const MOI = MathOptInterface @testset "MOI.FileFormats tests" begin - @testset "$(file)" for file in ["CBF", "LP", "MOF", "MPS", "SDPA"] + @testset "$(file)" for file in ["CBF", "LP", "MOF", "MPS", "NL", "SDPA"] include(joinpath(@__DIR__, file, "$(file).jl")) end @@ -66,6 +66,7 @@ const MOI = MathOptInterface (MOI.FileFormats.FORMAT_LP, MOI.FileFormats.LP.Model()), (MOI.FileFormats.FORMAT_MOF, MOI.FileFormats.MOF.Model()), (MOI.FileFormats.FORMAT_MPS, MOI.FileFormats.MPS.Model()), + (MOI.FileFormats.FORMAT_NL, MOI.FileFormats.NL.Model()), (MOI.FileFormats.FORMAT_SDPA, MOI.FileFormats.SDPA.Model()), ] @test typeof( @@ -83,6 +84,7 @@ const MOI = MathOptInterface (".lp", MOI.FileFormats.LP.Model()), (".mof.json", MOI.FileFormats.MOF.Model()), (".mps", MOI.FileFormats.MPS.Model()), + (".nl", MOI.FileFormats.NL.Model()), (".sdpa", MOI.FileFormats.SDPA.Model()), ] @test typeof(MOI.FileFormats.Model(filename = "a$(ext)")) == diff --git a/test/FileFormats/NL/NL.jl b/test/FileFormats/NL/NL.jl new file mode 100644 index 0000000000..a5f70182dc --- /dev/null +++ b/test/FileFormats/NL/NL.jl @@ -0,0 +1,759 @@ +module TestNLModel + +using MathOptInterface +const MOI = MathOptInterface +const NL = MOI.FileFormats.NL + +using Test + +function _test_nlexpr(expr::NL._NLExpr, nonlinear_terms, linear_terms, constant) + @test expr.is_linear == (length(nonlinear_terms) == 0) + @test expr.nonlinear_terms == nonlinear_terms + @test expr.linear_terms == linear_terms + @test expr.constant == constant + return +end +_test_nlexpr(x, args...) = _test_nlexpr(NL._NLExpr(x), args...) + +function test_nlexpr_singlevariable() + x = MOI.VariableIndex(1) + _test_nlexpr(MOI.SingleVariable(x), NL._NLTerm[], Dict(x => 1.0), 0.0) + return +end + +function test_nlexpr_scalaraffine() + x = MOI.VariableIndex.(1:3) + f = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.(1.0, x), 4.0) + return _test_nlexpr(f, NL._NLTerm[], Dict(x .=> 1), 4.0) +end + +function test_nlexpr_scalarquadratic() + x = MOI.VariableIndex(1) + f = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(1.1, x)], + [MOI.ScalarQuadraticTerm(2.0, x, x)], + 3.0, + ) + terms = [NL.OPMULT, x, x] + return _test_nlexpr(f, terms, Dict(x => 1.1), 3.0) +end + +function test_nlexpr_unary_addition() + x = MOI.VariableIndex(1) + return _test_nlexpr(:(+$x), [x], Dict(x => 0), 0.0) +end + +function test_nlexpr_binary_addition() + x = MOI.VariableIndex(1) + y = MOI.VariableIndex(2) + return _test_nlexpr( + :($x + $y), + [NL.OPPLUS, x, y], + Dict(x => 0.0, y => 0.0), + 0.0, + ) +end + +function test_nlexpr_nary_addition() + x = MOI.VariableIndex(1) + y = MOI.VariableIndex(2) + return _test_nlexpr( + :($x + $y + 1.0), + [NL.OPSUMLIST, 3, x, y, 1.0], + Dict(x => 0.0, y => 0.0), + 0.0, + ) +end + +function test_nlexpr_unary_subtraction() + x = MOI.VariableIndex(1) + return _test_nlexpr(:(-$x), [NL.OPUMINUS, x], Dict(x => 0.0), 0.0) +end + +function test_nlexpr_nary_multiplication() + x = MOI.VariableIndex(1) + return _test_nlexpr( + :($x * $x * 2.0), + [NL.OPMULT, x, NL.OPMULT, x, 2.0], + Dict(x => 0.0), + 0.0, + ) +end + +function test_nlexpr_unary_specialcase() + x = MOI.VariableIndex(1) + return _test_nlexpr( + :(cbrt($x)), + [NL.OPPOW, x, NL.OPDIV, 1, 3], + Dict(x => 0.0), + 0.0, + ) +end + +function test_nlexpr_binary_specialcase() + x = MOI.VariableIndex(1) + y = MOI.VariableIndex(2) + return _test_nlexpr( + :(\($x, $y)), + [NL.OPDIV, y, x], + Dict(x => 0.0, y => 0.0), + 0.0, + ) +end + +function test_nlexpr_unsupportedoperation() + x = MOI.VariableIndex(1) + err = ErrorException("Unsupported operation foo") + @test_throws err NL._NLExpr(:(foo($x))) + return +end + +function test_nlexpr_unsupportedexpression() + x = MOI.VariableIndex(1) + expr = :(1 <= $x <= 2) + err = ErrorException("Unsupported expression: $(expr)") + @test_throws err NL._NLExpr(expr) + return +end + +function test_nlexpr_ref() + x = MOI.VariableIndex(1) + return _test_nlexpr(:(x[$x]), [x], Dict(x => 0.0), 0.0) +end + +function test_nlconstraint_interval() + x = MOI.VariableIndex(1) + expr = :(1.0 <= $x <= 2.0) + con = NL._NLConstraint(expr, MOI.NLPBoundsPair(1.0, 2.0)) + @test con.lower == 1 + @test con.upper == 2 + @test con.opcode == 0 + @test con.expr == NL._NLExpr(expr.args[3]) +end + +function test_nlconstraint_lessthan() + x = MOI.VariableIndex(1) + expr = :($x <= 2.0) + con = NL._NLConstraint(expr, MOI.NLPBoundsPair(-Inf, 2.0)) + @test con.lower == -Inf + @test con.upper == 2 + @test con.opcode == 1 + @test con.expr == NL._NLExpr(expr.args[2]) +end + +function test_nlconstraint_greaterthan() + x = MOI.VariableIndex(1) + expr = :($x >= 2.0) + con = NL._NLConstraint(expr, MOI.NLPBoundsPair(2.0, Inf)) + @test con.lower == 2 + @test con.upper == Inf + @test con.opcode == 2 + @test con.expr == NL._NLExpr(expr.args[2]) +end + +function test_nlconstraint_equalto() + x = MOI.VariableIndex(1) + expr = :($x == 2.0) + con = NL._NLConstraint(expr, MOI.NLPBoundsPair(2.0, 2.0)) + @test con.lower == 2 + @test con.upper == 2 + @test con.opcode == 4 + @test con.expr == NL._NLExpr(expr.args[2]) +end + +function test_nlconstraint_interval_warn() + x = MOI.VariableIndex(1) + expr = :(2.0 <= $x <= 1.0) + @test_logs (:warn,) NL._NLConstraint(expr, MOI.NLPBoundsPair(1.0, 2.0)) +end + +function test_nlconstraint_lessthan_warn() + x = MOI.VariableIndex(1) + expr = :($x <= 1.0) + @test_logs (:warn,) NL._NLConstraint(expr, MOI.NLPBoundsPair(-Inf, 2.0)) +end + +function test_nlconstraint_greaterthan_warn() + x = MOI.VariableIndex(1) + expr = :($x >= 0.0) + @test_logs (:warn,) NL._NLConstraint(expr, MOI.NLPBoundsPair(1.0, Inf)) +end + +function test_nlconstraint_equalto_warn() + x = MOI.VariableIndex(1) + expr = :($x == 2.0) + @test_logs (:warn,) NL._NLConstraint(expr, MOI.NLPBoundsPair(1.0, 1.0)) +end + +function test_nlmodel_hs071() + model = MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()) + v = MOI.add_variables(model, 4) + l = [1.1, 1.2, 1.3, 1.4] + u = [5.1, 5.2, 5.3, 5.4] + start = [2.1, 2.2, 2.3, 2.4] + MOI.add_constraint.(model, MOI.SingleVariable.(v), MOI.GreaterThan.(l)) + MOI.add_constraint.(model, MOI.SingleVariable.(v), MOI.LessThan.(u)) + MOI.set.(model, MOI.VariablePrimalStart(), v, start) + lb, ub = [25.0, 40.0], [Inf, 40.0] + evaluator = MOI.Test.HS071(true) + block_data = MOI.NLPBlockData(MOI.NLPBoundsPair.(lb, ub), evaluator, true) + MOI.set(model, MOI.NLPBlock(), block_data) + MOI.set(model, MOI.ObjectiveSense(), MOI.MIN_SENSE) + n = NL.Model() + MOI.copy_to(n, model) + @test n.sense == MOI.MIN_SENSE + @test n.f == NL._NLExpr(MOI.objective_expr(evaluator)) + _test_nlexpr( + n.g[1].expr, + [NL.OPMULT, v[1], NL.OPMULT, v[2], NL.OPMULT, v[3], v[4]], + Dict(v .=> 0.0), + 0.0, + ) + @test n.g[1].lower == 25.0 + @test n.g[1].upper == Inf + @test n.g[1].opcode == 2 + _test_nlexpr( + n.g[2].expr, + [ + NL.OPSUMLIST, + 4, + NL.OPPOW, + v[1], + 2, + NL.OPPOW, + v[2], + 2, + NL.OPPOW, + v[3], + 2, + NL.OPPOW, + v[4], + 2, + ], + Dict(v .=> 0.0), + 0.0, + ) + @test n.g[2].lower == 40.0 + @test n.g[2].upper == 40.0 + @test n.g[2].opcode == 4 + @test length(n.h) == 0 + for i in 1:4 + @test n.x[v[i]].lower == l[i] + @test n.x[v[i]].upper == u[i] + @test n.x[v[i]].type == NL._CONTINUOUS + @test n.x[v[i]].jacobian_count == 2 + @test n.x[v[i]].in_nonlinear_constraint + @test n.x[v[i]].in_nonlinear_objective + @test 0 <= n.x[v[i]].order <= 3 + end + @test length(n.types[1]) == 4 + @test sprint(write, n) == """ + g3 1 1 0 + 4 2 1 0 1 0 + 2 1 + 0 0 + 4 4 4 + 0 0 0 1 + 0 0 0 0 0 + 8 4 + 0 0 + 0 0 0 0 0 + C0 + o2 + v0 + o2 + v1 + o2 + v2 + v3 + C1 + o54 + 4 + o5 + v0 + n2 + o5 + v1 + n2 + o5 + v2 + n2 + o5 + v3 + n2 + O0 0 + o0 + o2 + v0 + o2 + v3 + o54 + 3 + v0 + v1 + v2 + v2 + x4 + 0 2.1 + 1 2.2 + 2 2.3 + 3 2.4 + r + 2 25 + 4 40 + b + 0 1.1 5.1 + 0 1.2 5.2 + 0 1.3 5.3 + 0 1.4 5.4 + k3 + 2 + 4 + 6 + J0 4 + 0 0 + 1 0 + 2 0 + 3 0 + J1 4 + 0 0 + 1 0 + 2 0 + 3 0 + G0 4 + 0 0 + 1 0 + 2 0 + 3 0 + """ + return +end + +function test_nlmodel_hs071_linear_obj() + model = MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()) + v = MOI.add_variables(model, 4) + l = [1.1, 1.2, 1.3, 1.4] + u = [5.1, 5.2, 5.3, 5.4] + start = [2.1, 2.2, 2.3, 2.4] + MOI.add_constraint.(model, MOI.SingleVariable.(v), MOI.GreaterThan.(l)) + MOI.add_constraint.(model, MOI.SingleVariable.(v), MOI.LessThan.(u)) + MOI.add_constraint(model, MOI.SingleVariable(v[2]), MOI.ZeroOne()) + MOI.add_constraint(model, MOI.SingleVariable(v[3]), MOI.Integer()) + MOI.set.(model, MOI.VariablePrimalStart(), v, start) + lb, ub = [25.0, 40.0], [Inf, 40.0] + evaluator = MOI.Test.HS071(true) + block_data = MOI.NLPBlockData(MOI.NLPBoundsPair.(lb, ub), evaluator, false) + MOI.set(model, MOI.NLPBlock(), block_data) + f = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.(l, v), 2.0) + MOI.set(model, MOI.ObjectiveFunction{typeof(f)}(), f) + MOI.set(model, MOI.ObjectiveSense(), MOI.MAX_SENSE) + n = NL.Model() + MOI.copy_to(n, model) + @test n.sense == MOI.MAX_SENSE + @test n.f == NL._NLExpr(f) + _test_nlexpr( + n.g[1].expr, + [NL.OPMULT, v[1], NL.OPMULT, v[2], NL.OPMULT, v[3], v[4]], + Dict(v .=> 0.0), + 0.0, + ) + @test n.g[1].lower == 25.0 + @test n.g[1].upper == Inf + @test n.g[1].opcode == 2 + _test_nlexpr( + n.g[2].expr, + [ + NL.OPSUMLIST, + 4, + NL.OPPOW, + v[1], + 2, + NL.OPPOW, + v[2], + 2, + NL.OPPOW, + v[3], + 2, + NL.OPPOW, + v[4], + 2, + ], + Dict(v .=> 0.0), + 0.0, + ) + @test n.g[2].lower == 40.0 + @test n.g[2].upper == 40.0 + @test n.g[2].opcode == 4 + @test length(n.h) == 0 + types = [NL._CONTINUOUS, NL._BINARY, NL._INTEGER, NL._CONTINUOUS] + u[2] = 1.0 + for i in 1:4 + @test n.x[v[i]].lower == l[i] + @test n.x[v[i]].upper == u[i] + @test n.x[v[i]].type == types[i] + @test n.x[v[i]].jacobian_count == 2 + @test n.x[v[i]].in_nonlinear_constraint + @test !n.x[v[i]].in_nonlinear_objective + @test 0 <= n.x[v[i]].order <= 3 + end + @test length(n.types[3]) == 2 + @test length(n.types[4]) == 2 + @test v[1] in n.types[3] + @test v[2] in n.types[4] + @test v[3] in n.types[4] + @test v[4] in n.types[3] + @test sprint(write, n) == """ + g3 1 1 0 + 4 2 1 0 1 0 + 2 1 + 0 0 + 4 0 0 + 0 0 0 1 + 0 0 0 2 0 + 8 4 + 0 0 + 0 0 0 0 0 + C0 + o2 + v0 + o2 + v2 + o2 + v3 + v1 + C1 + o54 + 4 + o5 + v0 + n2 + o5 + v2 + n2 + o5 + v3 + n2 + o5 + v1 + n2 + O0 1 + n2 + x4 + 0 2.1 + 1 2.4 + 2 2.2 + 3 2.3 + r + 2 25 + 4 40 + b + 0 1.1 5.1 + 0 1.4 5.4 + 0 1.2 1 + 0 1.3 5.3 + k3 + 2 + 4 + 6 + J0 4 + 0 0 + 1 0 + 2 0 + 3 0 + J1 4 + 0 0 + 1 0 + 2 0 + 3 0 + G0 4 + 0 1.1 + 1 1.4 + 2 1.2 + 3 1.3 + """ + return +end + +function test_nlmodel_linear_quadratic() + model = MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()) + x = MOI.add_variables(model, 4) + MOI.add_constraint.(model, MOI.SingleVariable.(x), MOI.GreaterThan(0.0)) + MOI.add_constraint.(model, MOI.SingleVariable.(x), MOI.LessThan(2.0)) + MOI.add_constraint(model, MOI.SingleVariable(x[2]), MOI.ZeroOne()) + MOI.add_constraint(model, MOI.SingleVariable(x[3]), MOI.Integer()) + f = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.(1.0, x[2:4]), 2.0) + g = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(1.0, x[1])], + [MOI.ScalarQuadraticTerm(2.0, x[1], x[2])], + 3.0, + ) + h = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(1.0, x[3])], + [MOI.ScalarQuadraticTerm(1.0, x[1], x[2])], + 0.0, + ) + MOI.add_constraint(model, f, MOI.Interval(1.0, 10.0)) + MOI.add_constraint(model, g, MOI.LessThan(5.0)) + MOI.set(model, MOI.ObjectiveFunction{typeof(h)}(), h) + MOI.set(model, MOI.ObjectiveSense(), MOI.MAX_SENSE) + n = NL.Model() + MOI.copy_to(n, model) + @test n.sense == MOI.MAX_SENSE + @test n.f == NL._NLExpr(h) + terms = [NL.OPMULT, 2.0, NL.OPMULT, x[1], x[2]] + _test_nlexpr(n.g[1].expr, terms, Dict(x[1] => 1.0, x[2] => 0.0), 3.0) + @test n.g[1].opcode == 1 + @test n.g[1].lower == -Inf + @test n.g[1].upper == 5.0 + @test n.h[1].expr == NL._NLExpr(f) + @test n.h[1].opcode == 0 + @test n.h[1].lower == 1.0 + @test n.h[1].upper == 10.0 + @test n.types[1] == [x[1]] # Continuous in both + @test n.types[2] == [x[2]] # Discrete in both + @test n.types[6] == [x[3]] # Discrete in objective only + @test n.types[7] == [x[4]] # Continuous in linear + @test sprint(write, n) == """ + g3 1 1 0 + 4 2 1 1 0 0 + 1 1 + 0 0 + 2 3 2 + 0 0 0 1 + 0 0 1 0 1 + 5 3 + 0 0 + 0 0 0 0 0 + C0 + o0 + n3 + o2 + n2 + o2 + v0 + v1 + C1 + n2 + O0 1 + o2 + v0 + v1 + x4 + 0 0 + 1 0 + 2 0 + 3 0 + r + 1 5 + 0 1 10 + b + 0 0 2 + 0 0 1 + 0 0 2 + 0 0 2 + k3 + 1 + 3 + 4 + J0 2 + 0 1 + 1 0 + J1 3 + 1 1 + 2 1 + 3 1 + G0 3 + 0 0 + 1 0 + 2 1 + """ + return +end + +function test_nlmodel_quadratic_interval() + model = MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()) + x = MOI.add_variable(model) + g = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(1.0, x)], + [MOI.ScalarQuadraticTerm(2.0, x, x)], + 3.0, + ) + MOI.add_constraint(model, g, MOI.Interval(1.0, 10.0)) + n = NL.Model() + MOI.copy_to(n, model) + @test sprint(write, n) == """ + g3 1 1 0 + 1 1 1 1 0 0 + 1 1 + 0 0 + 1 0 0 + 0 0 0 1 + 0 0 0 0 0 + 1 0 + 0 0 + 0 0 0 0 0 + C0 + o0 + n3 + o2 + v0 + v0 + O0 0 + n0 + x1 + 0 0 + r + 0 1 10 + b + 3 + k0 + J0 1 + 0 1 + """ + return +end + +function test_eval_singlevariable() + x = MOI.VariableIndex(1) + f = NL._NLExpr(MOI.SingleVariable(x)) + @test NL._evaluate(f, Dict(x => 1.2)) == 1.2 +end + +function test_eval_scalaraffine() + x = MOI.VariableIndex.(1:3) + f = NL._NLExpr(MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.(1.0, x), 4.0)) + @test NL._evaluate(f, Dict(x[i] => Float64(i) for i in 1:3)) == 10.0 +end + +function test_eval_scalarquadratic() + x = MOI.VariableIndex(1) + f = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(1.1, x)], + [MOI.ScalarQuadraticTerm(2.0, x, x)], + 3.0, + ) + @test NL._evaluate(NL._NLExpr(f), Dict(x => 1.1)) == 5.42 +end + +function test_eval_unary_addition() + x = MOI.VariableIndex(1) + @test NL._evaluate(NL._NLExpr(:(+$x)), Dict(x => 1.1)) == 1.1 +end + +function test_eval_binary_addition() + x = MOI.VariableIndex(1) + y = MOI.VariableIndex(2) + @test NL._evaluate(NL._NLExpr(:($x + $y)), Dict(x => 1.1, y => 2.2)) ≈ 3.3 +end + +function test_eval_nary_addition() + x = MOI.VariableIndex(1) + y = MOI.VariableIndex(2) + @test NL._evaluate(NL._NLExpr(:($x + $y + 1.0)), Dict(x => 1.1, y => 2.2)) ≈ + 4.3 +end + +function test_eval_unary_subtraction() + x = MOI.VariableIndex(1) + @test NL._evaluate(NL._NLExpr(:(-$x)), Dict(x => 1.1)) == -1.1 +end + +function test_eval_nary_multiplication() + x = MOI.VariableIndex(1) + @test NL._evaluate(NL._NLExpr(:($x * $x * 2.0)), Dict(x => 1.1)) ≈ 2.42 +end + +function test_eval_unary_specialcases() + x = MOI.VariableIndex(1) + S = [:acoth, :asec, :acsc, :acot, :asecd, :acscd, :acotd] + for (k, v) in NL._UNARY_SPECIAL_CASES + expr = NL._NLExpr(:($k($x))) + xv = k in S ? 1.49 : 0.49 + @test NL._evaluate(expr, Dict(x => xv)) ≈ getfield(Main, k)(xv) + end +end + +function test_eval_binary_specialcases() + x = MOI.VariableIndex(1) + y = MOI.VariableIndex(2) + for (k, v) in NL._BINARY_SPECIAL_CASES + expr = NL._NLExpr(:($k($x, $y))) + target = getfield(Main, k)(1.0, 2.0) + @test NL._evaluate(expr, Dict(x => 1.0, y => 2.0)) ≈ target + end +end + +""" + test_issue_79() + +Test the problem + + min (z - 0.5)^2 = z^2 - z + 1/4 + s.t. x * z <= 0 + z ∈ {0, 1} +""" +function test_issue_79() + model = MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()) + x = MOI.add_variable(model) + z = MOI.add_variable(model) + MOI.add_constraint(model, MOI.SingleVariable(z), MOI.ZeroOne()) + f = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(-1.0, z)], + [MOI.ScalarQuadraticTerm(1.0, z, z)], + 0.25, + ) + MOI.set(model, MOI.ObjectiveFunction{typeof(f)}(), f) + MOI.set(model, MOI.ObjectiveSense(), MOI.MIN_SENSE) + MOI.add_constraint( + model, + MOI.ScalarQuadraticFunction( + MOI.ScalarAffineTerm{Float64}[], + [MOI.ScalarQuadraticTerm(1.0, x, z)], + 0.0, + ), + MOI.LessThan(0.0), + ) + n = NL.Model() + MOI.copy_to(n, model) + # x is continuous in a nonlinear constraint [Priority 3] + # z is discrete in a nonlinear constraint and objective [Priority 2] + @test n.x[x].order == 1 + @test n.x[z].order == 0 +end + +function test_malformed_constraint_error() + model = MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()) + x = MOI.add_variable(model) + MOI.add_constraint( + model, + MOI.ScalarAffineFunction(MOI.ScalarAffineTerm{Float64}[], 1.0), + MOI.LessThan(0.0), + ) + n = NL.Model() + @test_throws ErrorException MOI.copy_to(n, model) +end + +struct NoExprGraph <: MOI.AbstractNLPEvaluator end +MOI.features_available(::NoExprGraph) = Symbol[] + +function test_noexpr_graph() + model = MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()) + block_data = MOI.NLPBlockData(MOI.NLPBoundsPair[], NoExprGraph(), false) + MOI.set(model, MOI.NLPBlock(), block_data) + n = NL.Model() + @test_throws ErrorException MOI.copy_to(n, model) +end + +function runtests() + for name in names(@__MODULE__; all = true) + if startswith("$(name)", "test_") + @testset "$(name)" begin + getfield(@__MODULE__, name)() + end + end + end +end + +end + +TestNLModel.runtests() From e1224f0471882b2f74466a502e10b2bf48d01df9 Mon Sep 17 00:00:00 2001 From: odow Date: Fri, 30 Apr 2021 16:55:52 +1200 Subject: [PATCH 02/12] Add gen.jl --- src/FileFormats/NL/gen.jl | 15 +++++++++++++++ 1 file changed, 15 insertions(+) create mode 100644 src/FileFormats/NL/gen.jl diff --git a/src/FileFormats/NL/gen.jl b/src/FileFormats/NL/gen.jl new file mode 100644 index 0000000000..c14365929a --- /dev/null +++ b/src/FileFormats/NL/gen.jl @@ -0,0 +1,15 @@ +# This script builds the list of recognized ASL opcodes using the header files +# in ASL_jll. Only re-run it if AMPL adds new opcodes (which is unlikely). + +using ASL_jll + +open("opcode.jl", "w") do io + println(""" + # Do not modify. This file is automatically created by the script in `gen.jl`. + """) + filename = joinpath(ASL_jll.artifact_dir, "include", "opcode.hd") + for line in readlines(filename) + items = split(line, c -> c == '\t' || c == ' '; keepempty = false) + println(io, "const ", items[2], " = ", items[3]) + end +end From 0b88c2b097ce36b1eba92866bb28ad07f6ed1fcc Mon Sep 17 00:00:00 2001 From: Oscar Dowson Date: Fri, 30 Apr 2021 17:00:21 +1200 Subject: [PATCH 03/12] Update gen.jl --- src/FileFormats/NL/gen.jl | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/src/FileFormats/NL/gen.jl b/src/FileFormats/NL/gen.jl index c14365929a..3c440e975f 100644 --- a/src/FileFormats/NL/gen.jl +++ b/src/FileFormats/NL/gen.jl @@ -4,9 +4,10 @@ using ASL_jll open("opcode.jl", "w") do io - println(""" - # Do not modify. This file is automatically created by the script in `gen.jl`. - """) + println( + "# Do not modify. This file is automatically created by the script " * + "in `gen.jl`.\n", + ) filename = joinpath(ASL_jll.artifact_dir, "include", "opcode.hd") for line in readlines(filename) items = split(line, c -> c == '\t' || c == ' '; keepempty = false) From 6bcedd54f406dc57896bb7b3ef2cb19d813eaa3c Mon Sep 17 00:00:00 2001 From: odow Date: Fri, 7 May 2021 13:29:07 +1200 Subject: [PATCH 04/12] Updates and better tests --- src/FileFormats/NL/NL.jl | 9 ++--- src/FileFormats/NL/NLExpr.jl | 14 ++++---- src/FileFormats/NL/gen.jl | 6 +++- test/FileFormats/NL/NL.jl | 65 ++++++++++++++++++++++++++++++++++++ 4 files changed, 82 insertions(+), 12 deletions(-) diff --git a/src/FileFormats/NL/NL.jl b/src/FileFormats/NL/NL.jl index 0c77177adf..8ae1242480 100644 --- a/src/FileFormats/NL/NL.jl +++ b/src/FileFormats/NL/NL.jl @@ -259,8 +259,6 @@ struct _LinearNLPEvaluator <: MOI.AbstractNLPEvaluator end MOI.features_available(::_LinearNLPEvaluator) = [:ExprGraph] MOI.initialize(::_LinearNLPEvaluator, ::Vector{Symbol}) = nothing -MOI.Utilities.supports_default_copy_to(::Model, ::Bool) = false - function MOI.copy_to( dest::Model, model::MOI.ModelLike; @@ -296,10 +294,13 @@ function MOI.copy_to( ) end dest.nlpblock_dim = length(dest.g) + starts = MOI.supports(model, MOI.VariablePrimalStart(), MOI.VariableIndex) for x in MOI.get(model, MOI.ListOfVariableIndices()) dest.x[x] = _VariableInfo() - start = MOI.get(model, MOI.VariablePrimalStart(), x) - MOI.set(dest, MOI.VariablePrimalStart(), x, start) + if starts + start = MOI.get(model, MOI.VariablePrimalStart(), x) + MOI.set(dest, MOI.VariablePrimalStart(), x, start) + end mapping[x] = x end dest.sense = MOI.get(model, MOI.ObjectiveSense()) diff --git a/src/FileFormats/NL/NLExpr.jl b/src/FileFormats/NL/NLExpr.jl index 5dc39e3371..9fe1b95249 100644 --- a/src/FileFormats/NL/NLExpr.jl +++ b/src/FileFormats/NL/NLExpr.jl @@ -310,7 +310,7 @@ end function _NLExpr(x::MOI.ScalarAffineFunction) linear = Dict{MOI.VariableIndex,Float64}() for (i, term) in enumerate(x.terms) - _add_or_set(linear, term.variable_index, term.coefficient) + _add_or_set(linear, term.variable, term.coefficient) end return _NLExpr(true, _NLTerm[], linear, x.constant) end @@ -318,7 +318,7 @@ end function _NLExpr(x::MOI.ScalarQuadraticFunction) linear = Dict{MOI.VariableIndex,Float64}() for (i, term) in enumerate(x.affine_terms) - _add_or_set(linear, term.variable_index, term.coefficient) + _add_or_set(linear, term.variable, term.coefficient) end terms = _NLTerm[] N = length(x.quadratic_terms) @@ -335,7 +335,7 @@ function _NLExpr(x::MOI.ScalarQuadraticFunction) for term in x.quadratic_terms coefficient = term.coefficient # MOI defines quadratic as 1/2 x' Q x :( - if term.variable_index_1 == term.variable_index_2 + if term.variable_1 == term.variable_2 coefficient *= 0.5 end # Optimization: no need for the OPMULT if the coefficient is 1. @@ -344,12 +344,12 @@ function _NLExpr(x::MOI.ScalarQuadraticFunction) push!(terms, coefficient) end push!(terms, OPMULT) - push!(terms, term.variable_index_1) - push!(terms, term.variable_index_2) + push!(terms, term.variable_1) + push!(terms, term.variable_2) # For the Jacobian sparsity patterns, we need to add a linear term, even # if the variable only appears nonlinearly. - _add_or_set(linear, term.variable_index_1, 0.0) - _add_or_set(linear, term.variable_index_2, 0.0) + _add_or_set(linear, term.variable_1, 0.0) + _add_or_set(linear, term.variable_2, 0.0) end return _NLExpr(false, terms, linear, x.constant) end diff --git a/src/FileFormats/NL/gen.jl b/src/FileFormats/NL/gen.jl index 3c440e975f..ed3321b864 100644 --- a/src/FileFormats/NL/gen.jl +++ b/src/FileFormats/NL/gen.jl @@ -1,5 +1,9 @@ # This script builds the list of recognized ASL opcodes using the header files -# in ASL_jll. Only re-run it if AMPL adds new opcodes (which is unlikely). +# in ASL_jll. +# +# Because it is used for code generation, it is not included by +# `MathOptInterface.FileFormats.NL`. Only re-run it if AMPL adds new opcodes +# (which is unlikely). using ASL_jll diff --git a/test/FileFormats/NL/NL.jl b/test/FileFormats/NL/NL.jl index a5f70182dc..0f6978ddfd 100644 --- a/test/FileFormats/NL/NL.jl +++ b/test/FileFormats/NL/NL.jl @@ -197,10 +197,14 @@ function test_nlmodel_hs071() lb, ub = [25.0, 40.0], [Inf, 40.0] evaluator = MOI.Test.HS071(true) block_data = MOI.NLPBlockData(MOI.NLPBoundsPair.(lb, ub), evaluator, true) + @test MOI.supports(model, MOI.NLPBlock()) MOI.set(model, MOI.NLPBlock(), block_data) + @test MOI.supports(model, MOI.ObjectiveSense()) MOI.set(model, MOI.ObjectiveSense(), MOI.MIN_SENSE) n = NL.Model() + @test MOI.is_empty(n) MOI.copy_to(n, model) + @test !MOI.is_empty(n) @test n.sense == MOI.MIN_SENSE @test n.f == NL._NLExpr(MOI.objective_expr(evaluator)) _test_nlexpr( @@ -339,12 +343,14 @@ function test_nlmodel_hs071_linear_obj() MOI.add_constraint.(model, MOI.SingleVariable.(v), MOI.LessThan.(u)) MOI.add_constraint(model, MOI.SingleVariable(v[2]), MOI.ZeroOne()) MOI.add_constraint(model, MOI.SingleVariable(v[3]), MOI.Integer()) + @test MOI.supports(model, MOI.VariablePrimalStart(), MOI.VariableIndex) MOI.set.(model, MOI.VariablePrimalStart(), v, start) lb, ub = [25.0, 40.0], [Inf, 40.0] evaluator = MOI.Test.HS071(true) block_data = MOI.NLPBlockData(MOI.NLPBoundsPair.(lb, ub), evaluator, false) MOI.set(model, MOI.NLPBlock(), block_data) f = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.(l, v), 2.0) + @test MOI.supports(model, MOI.ObjectiveFunction{typeof(f)}()) MOI.set(model, MOI.ObjectiveFunction{typeof(f)}(), f) MOI.set(model, MOI.ObjectiveSense(), MOI.MAX_SENSE) n = NL.Model() @@ -744,6 +750,65 @@ function test_noexpr_graph() @test_throws ErrorException MOI.copy_to(n, model) end +function test_Name() + model = NL.Model() + @test MOI.supports(model, MOI.Name()) + MOI.set(model, MOI.Name(), "MyModel") + @test MOI.get(model, MOI.Name()) == "MyModel" +end + +function test_RawParameter() + model = NL.Model() + @test MOI.supports(model, MOI.RawParameter("print_level")) + MOI.set(model, MOI.RawParameter("print_level"), 0) + @test MOI.get(model, MOI.RawParameter("print_level")) == 0 +end + +function test_SolverName() + model = NL.Model() + @test MOI.get(model, MOI.SolverName()) == "AmplNLWriter" +end + +function test_show() + model = NL.Model() + @test sprint(show, model) == "An AMPL (.nl) model" +end + +function test_linear_constraint_types() + model = MOI.Utilities.Model{Float64}() + y = MOI.add_variables(model, 3) + MOI.add_constraint(model, MOI.SingleVariable(y[1]), MOI.ZeroOne()) + MOI.add_constraint(model, MOI.SingleVariable(y[2]), MOI.Integer()) + for set in [ + MOI.GreaterThan(0.0), + MOI.LessThan(1.0), + MOI.EqualTo(2.0), + MOI.Interval(3.0, 4.0), + ] + x = MOI.add_variable(model) + MOI.add_constraint(model, MOI.SingleVariable(x), set) + MOI.add_constraint( + model, + MOI.ScalarAffineFunction([MOI.ScalarAffineTerm(1.0, x)], 0.0), + set, + ) + end + n = NL.Model() + map = MOI.copy_to(n, model) + @test length(n.g) == 0 + @test length(n.h) == 4 + @test length(n.x) == 7 + @test n.order == [ + MOI.VariableIndex(3), + MOI.VariableIndex(4), + MOI.VariableIndex(5), + MOI.VariableIndex(6), + MOI.VariableIndex(7), + MOI.VariableIndex(1), + MOI.VariableIndex(2), + ] +end + function runtests() for name in names(@__MODULE__; all = true) if startswith("$(name)", "test_") From 97e611962407c3315697d1a4922bf6c32bce6f48 Mon Sep 17 00:00:00 2001 From: odow Date: Fri, 7 May 2021 15:04:09 +1200 Subject: [PATCH 05/12] Fix tests --- src/FileFormats/NL/NL.jl | 11 +---------- test/FileFormats/NL/NL.jl | 7 ------- 2 files changed, 1 insertion(+), 17 deletions(-) diff --git a/src/FileFormats/NL/NL.jl b/src/FileFormats/NL/NL.jl index 8ae1242480..56ddeeb907 100644 --- a/src/FileFormats/NL/NL.jl +++ b/src/FileFormats/NL/NL.jl @@ -220,15 +220,6 @@ end MOI.supports(::Model, ::MOI.ObjectiveSense) = true MOI.supports(::Model, ::MOI.ObjectiveFunction{<:_SCALAR_FUNCTIONS}) = true -MOI.supports(::Model, ::MOI.RawParameter) = true -function MOI.get(model::Model, attr::MOI.RawParameter) - return get(model.options, attr.name, nothing) -end -function MOI.set(model::Model, attr::MOI.RawParameter, value) - model.options[attr.name] = value - return -end - # ============================================================================== function MOI.supports( @@ -306,7 +297,7 @@ function MOI.copy_to( dest.sense = MOI.get(model, MOI.ObjectiveSense()) resize!(dest.order, length(dest.x)) # Now deal with the normal MOI constraints. - for (F, S) in MOI.get(model, MOI.ListOfConstraints()) + for (F, S) in MOI.get(model, MOI.ListOfConstraintTypesPresents()) _process_constraint(dest, model, F, S, mapping) end # Correct bounds of binary variables. Mainly because AMPL doesn't have the diff --git a/test/FileFormats/NL/NL.jl b/test/FileFormats/NL/NL.jl index 0f6978ddfd..85a8c5e51a 100644 --- a/test/FileFormats/NL/NL.jl +++ b/test/FileFormats/NL/NL.jl @@ -757,13 +757,6 @@ function test_Name() @test MOI.get(model, MOI.Name()) == "MyModel" end -function test_RawParameter() - model = NL.Model() - @test MOI.supports(model, MOI.RawParameter("print_level")) - MOI.set(model, MOI.RawParameter("print_level"), 0) - @test MOI.get(model, MOI.RawParameter("print_level")) == 0 -end - function test_SolverName() model = NL.Model() @test MOI.get(model, MOI.SolverName()) == "AmplNLWriter" From 5a7e0f9d68abbcdd953b91a1b08e9596d88098b5 Mon Sep 17 00:00:00 2001 From: odow Date: Fri, 7 May 2021 16:13:38 +1200 Subject: [PATCH 06/12] Fix typo --- src/FileFormats/NL/NL.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/FileFormats/NL/NL.jl b/src/FileFormats/NL/NL.jl index 56ddeeb907..2a780616be 100644 --- a/src/FileFormats/NL/NL.jl +++ b/src/FileFormats/NL/NL.jl @@ -297,7 +297,7 @@ function MOI.copy_to( dest.sense = MOI.get(model, MOI.ObjectiveSense()) resize!(dest.order, length(dest.x)) # Now deal with the normal MOI constraints. - for (F, S) in MOI.get(model, MOI.ListOfConstraintTypesPresents()) + for (F, S) in MOI.get(model, MOI.ListOfConstraintTypesPresent()) _process_constraint(dest, model, F, S, mapping) end # Correct bounds of binary variables. Mainly because AMPL doesn't have the From f3582a2a000a2df023ab35eb40dc8aaeda3e250c Mon Sep 17 00:00:00 2001 From: odow Date: Fri, 7 May 2021 17:15:09 +1200 Subject: [PATCH 07/12] Improve coverage --- test/FileFormats/NL/NL.jl | 34 +++++++++++++++++++++++++++++----- 1 file changed, 29 insertions(+), 5 deletions(-) diff --git a/test/FileFormats/NL/NL.jl b/test/FileFormats/NL/NL.jl index 85a8c5e51a..a962c1b847 100644 --- a/test/FileFormats/NL/NL.jl +++ b/test/FileFormats/NL/NL.jl @@ -197,11 +197,11 @@ function test_nlmodel_hs071() lb, ub = [25.0, 40.0], [Inf, 40.0] evaluator = MOI.Test.HS071(true) block_data = MOI.NLPBlockData(MOI.NLPBoundsPair.(lb, ub), evaluator, true) - @test MOI.supports(model, MOI.NLPBlock()) MOI.set(model, MOI.NLPBlock(), block_data) - @test MOI.supports(model, MOI.ObjectiveSense()) MOI.set(model, MOI.ObjectiveSense(), MOI.MIN_SENSE) n = NL.Model() + @test MOI.supports(n, MOI.NLPBlock()) + @test MOI.supports(n, MOI.ObjectiveSense()) @test MOI.is_empty(n) MOI.copy_to(n, model) @test !MOI.is_empty(n) @@ -343,17 +343,17 @@ function test_nlmodel_hs071_linear_obj() MOI.add_constraint.(model, MOI.SingleVariable.(v), MOI.LessThan.(u)) MOI.add_constraint(model, MOI.SingleVariable(v[2]), MOI.ZeroOne()) MOI.add_constraint(model, MOI.SingleVariable(v[3]), MOI.Integer()) - @test MOI.supports(model, MOI.VariablePrimalStart(), MOI.VariableIndex) MOI.set.(model, MOI.VariablePrimalStart(), v, start) lb, ub = [25.0, 40.0], [Inf, 40.0] evaluator = MOI.Test.HS071(true) block_data = MOI.NLPBlockData(MOI.NLPBoundsPair.(lb, ub), evaluator, false) MOI.set(model, MOI.NLPBlock(), block_data) f = MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.(l, v), 2.0) - @test MOI.supports(model, MOI.ObjectiveFunction{typeof(f)}()) MOI.set(model, MOI.ObjectiveFunction{typeof(f)}(), f) MOI.set(model, MOI.ObjectiveSense(), MOI.MAX_SENSE) n = NL.Model() + @test MOI.supports(n, MOI.VariablePrimalStart(), MOI.VariableIndex) + @test MOI.supports(n, MOI.ObjectiveFunction{typeof(f)}()) MOI.copy_to(n, model) @test n.sense == MOI.MAX_SENSE @test n.f == NL._NLExpr(f) @@ -769,9 +769,12 @@ end function test_linear_constraint_types() model = MOI.Utilities.Model{Float64}() + n = NL.Model() y = MOI.add_variables(model, 3) MOI.add_constraint(model, MOI.SingleVariable(y[1]), MOI.ZeroOne()) MOI.add_constraint(model, MOI.SingleVariable(y[2]), MOI.Integer()) + @test MOI.supports_constraint(n, MOI.SingleVariable, MOI.ZeroOne) + @test MOI.supports_constraint(n, MOI.SingleVariable, MOI.Integer) for set in [ MOI.GreaterThan(0.0), MOI.LessThan(1.0), @@ -779,6 +782,16 @@ function test_linear_constraint_types() MOI.Interval(3.0, 4.0), ] x = MOI.add_variable(model) + @test MOI.supports_constraint( + n, + MOI.SingleVariable, + typeof(set), + ) + @test MOI.supports_constraint( + n, + MOI.ScalarAffineFunction{Float64}, + typeof(set), + ) MOI.add_constraint(model, MOI.SingleVariable(x), set) MOI.add_constraint( model, @@ -786,7 +799,6 @@ function test_linear_constraint_types() set, ) end - n = NL.Model() map = MOI.copy_to(n, model) @test length(n.g) == 0 @test length(n.h) == 4 @@ -800,6 +812,18 @@ function test_linear_constraint_types() MOI.VariableIndex(1), MOI.VariableIndex(2), ] + file = sprint(Base.write, n) +end + + +function test_empty() + model = MOI.Utilities.Model{Float64}() + y = MOI.add_variables(model, 3) + n = NL.Model() + map = MOI.copy_to(n, model) + @test !MOI.is_empty(model) + MOI.empty!(model) + @test MOI.is_empty(model) end function runtests() From fa79e5b6092d346bd72f1b9ae6dbc02741b0e2aa Mon Sep 17 00:00:00 2001 From: odow Date: Mon, 10 May 2021 08:40:59 +1200 Subject: [PATCH 08/12] Various fixes --- src/FileFormats/NL/NL.jl | 24 +++++++------- test/FileFormats/NL/NL.jl | 67 +++++++++++++++++++++++++++++++++++---- 2 files changed, 72 insertions(+), 19 deletions(-) diff --git a/src/FileFormats/NL/NL.jl b/src/FileFormats/NL/NL.jl index 2a780616be..1da9cba24e 100644 --- a/src/FileFormats/NL/NL.jl +++ b/src/FileFormats/NL/NL.jl @@ -103,7 +103,7 @@ mutable struct _VariableInfo jacobian_count::Int # If the variable appears in the objective. in_nonlinear_objective::Bool - # If the objetive appears in a nonlinear constraint. + # If the variable appears in a nonlinear constraint. in_nonlinear_constraint::Bool # The 0-indexed column of the variable. Computed right at the end. order::Int @@ -298,7 +298,7 @@ function MOI.copy_to( resize!(dest.order, length(dest.x)) # Now deal with the normal MOI constraints. for (F, S) in MOI.get(model, MOI.ListOfConstraintTypesPresent()) - _process_constraint(dest, model, F, S, mapping) + _process_constraint(dest, model, mapping.conmap[F, S]) end # Correct bounds of binary variables. Mainly because AMPL doesn't have the # concept of binary nonlinear variables, but it does have binary linear @@ -444,7 +444,11 @@ _set_to_bounds(set::MOI.LessThan) = (1, -Inf, set.upper) _set_to_bounds(set::MOI.GreaterThan) = (2, set.lower, Inf) _set_to_bounds(set::MOI.EqualTo) = (4, set.value, set.value) -function _process_constraint(dest::Model, model, F, S, mapping) +function _process_constraint( + dest::Model, + model, + mapping::MOI.Utilities.DoubleDicts.IndexWithType{F,S}, +) where {F,S} for ci in MOI.get(model, MOI.ListOfConstraintIndices{F,S}()) f = MOI.get(model, MOI.ConstraintFunction(), ci) s = MOI.get(model, MOI.ConstraintSet(), ci) @@ -456,7 +460,7 @@ function _process_constraint(dest::Model, model, F, S, mapping) else error( "Malformed constraint. There are no variables and the " * - "bounds don't make sense.", + "function constant $(con.expr.constant) is not in [$l, $u]", ) end elseif con.expr.is_linear @@ -473,10 +477,8 @@ end function _process_constraint( dest::Model, model, - F::Type{MOI.SingleVariable}, - S, - mapping, -) + mapping::MOI.Utilities.DoubleDicts.IndexWithType{F,S}, +) where {F<:MOI.SingleVariable,S} for ci in MOI.get(model, MOI.ListOfConstraintIndices{F,S}()) mapping[ci] = ci f = MOI.get(model, MOI.ConstraintFunction(), ci) @@ -495,10 +497,8 @@ end function _process_constraint( dest::Model, model, - F::Type{MOI.SingleVariable}, - S::Union{Type{MOI.ZeroOne},Type{MOI.Integer}}, - mapping, -) + mapping::MOI.Utilities.DoubleDicts.IndexWithType{F,S}, +) where {F<:MOI.SingleVariable,S<:Union{MOI.ZeroOne,MOI.Integer}} for ci in MOI.get(model, MOI.ListOfConstraintIndices{F,S}()) mapping[ci] = ci f = MOI.get(model, MOI.ConstraintFunction(), ci) diff --git a/test/FileFormats/NL/NL.jl b/test/FileFormats/NL/NL.jl index a962c1b847..87c00b733b 100644 --- a/test/FileFormats/NL/NL.jl +++ b/test/FileFormats/NL/NL.jl @@ -782,11 +782,7 @@ function test_linear_constraint_types() MOI.Interval(3.0, 4.0), ] x = MOI.add_variable(model) - @test MOI.supports_constraint( - n, - MOI.SingleVariable, - typeof(set), - ) + @test MOI.supports_constraint(n, MOI.SingleVariable, typeof(set)) @test MOI.supports_constraint( n, MOI.ScalarAffineFunction{Float64}, @@ -812,10 +808,67 @@ function test_linear_constraint_types() MOI.VariableIndex(1), MOI.VariableIndex(2), ] - file = sprint(Base.write, n) + @test sprint(Base.write, n) == """ + g3 1 1 0 + 7 4 1 1 1 0 + 0 1 + 0 0 + 0 0 0 + 0 0 0 1 + 1 1 0 0 0 + 4 0 + 0 0 + 0 0 0 0 0 + C0 + n0 + C1 + n0 + C2 + n0 + C3 + n0 + O0 0 + n0 + x7 + 0 0 + 1 0 + 2 0 + 3 0 + 4 0 + 5 0 + 6 0 + r + 4 2 + 2 0 + 1 1 + 0 3 4 + b + 3 + 2 0 + 1 1 + 4 2 + 0 3 4 + 0 0 1 + 3 + k6 + 0 + 1 + 2 + 3 + 4 + 4 + J0 1 + 3 1 + J1 1 + 1 1 + J2 1 + 2 1 + J3 1 + 4 1 + """ + return end - function test_empty() model = MOI.Utilities.Model{Float64}() y = MOI.add_variables(model, 3) From df326b1f3f51f74219ac6095ff7fbfcfd91a930e Mon Sep 17 00:00:00 2001 From: odow Date: Mon, 10 May 2021 11:03:29 +1200 Subject: [PATCH 09/12] Improve coverage --- src/FileFormats/NL/NL.jl | 10 --- src/FileFormats/NL/NLExpr.jl | 4 +- test/FileFormats/NL/NL.jl | 124 +++++++++++++++++++++++++++++++++-- 3 files changed, 118 insertions(+), 20 deletions(-) diff --git a/src/FileFormats/NL/NL.jl b/src/FileFormats/NL/NL.jl index 1da9cba24e..a5fe132518 100644 --- a/src/FileFormats/NL/NL.jl +++ b/src/FileFormats/NL/NL.jl @@ -162,16 +162,6 @@ MOI.get(model::Model, ::MOI.Name) = model.name MOI.set(model::Model, ::MOI.Name, name::String) = (model.name = name) function MOI.empty!(model::Model) - model.results = _NLResults( - "Optimize not called.", - MOI.OPTIMIZE_NOT_CALLED, - MOI.NO_SOLUTION, - NaN, - Dict{MOI.VariableIndex,Float64}(), - Float64[], - Dict{MOI.VariableIndex,Float64}(), - Dict{MOI.VariableIndex,Float64}(), - ) model.f = _NLExpr(false, _NLTerm[], Dict{MOI.VariableIndex,Float64}(), 0.0) empty!(model.g) model.nlpblock_dim = 0 diff --git a/src/FileFormats/NL/NLExpr.jl b/src/FileFormats/NL/NLExpr.jl index 9fe1b95249..e8ce2e58d3 100644 --- a/src/FileFormats/NL/NLExpr.jl +++ b/src/FileFormats/NL/NLExpr.jl @@ -393,10 +393,8 @@ function _process_expr!(expr::_NLExpr, arg::Expr) return _process_expr!(expr, arg.args[2]) elseif arg == :() return # Some evalators return a null objective of `:()`. - else - error("Unsupported expression: $(arg)") end - return + return error("Unsupported expression: $(arg)") end function _process_expr!(expr::_NLExpr, args::Vector{Any}) diff --git a/test/FileFormats/NL/NL.jl b/test/FileFormats/NL/NL.jl index 87c00b733b..ef0decb596 100644 --- a/test/FileFormats/NL/NL.jl +++ b/test/FileFormats/NL/NL.jl @@ -6,14 +6,27 @@ const NL = MOI.FileFormats.NL using Test -function _test_nlexpr(expr::NL._NLExpr, nonlinear_terms, linear_terms, constant) - @test expr.is_linear == (length(nonlinear_terms) == 0) +function _test_nlexpr( + expr::NL._NLExpr, + nonlinear_terms, + linear_terms, + constant; + force_nonlinear::Bool = false, +) + if force_nonlinear + @test expr.is_linear == false + else + @test expr.is_linear == (length(nonlinear_terms) == 0) + end @test expr.nonlinear_terms == nonlinear_terms @test expr.linear_terms == linear_terms @test expr.constant == constant return end -_test_nlexpr(x, args...) = _test_nlexpr(NL._NLExpr(x), args...) + +function _test_nlexpr(x, args...; kwargs...) + return _test_nlexpr(NL._NLExpr(x), args...; kwargs...) +end function test_nlexpr_singlevariable() x = MOI.VariableIndex(1) @@ -27,7 +40,23 @@ function test_nlexpr_scalaraffine() return _test_nlexpr(f, NL._NLTerm[], Dict(x .=> 1), 4.0) end -function test_nlexpr_scalarquadratic() +function test_nlexpr_scalarquadratic_0() + x = MOI.VariableIndex(1) + f = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(1.1, x)], + MOI.ScalarQuadraticTerm{Float64}[], + 3.0, + ) + return _test_nlexpr( + f, + NL._NLTerm[], + Dict(x => 1.1), + 3.0; + force_nonlinear = true, + ) +end + +function test_nlexpr_scalarquadratic_1a() x = MOI.VariableIndex(1) f = MOI.ScalarQuadraticFunction( [MOI.ScalarAffineTerm(1.1, x)], @@ -38,6 +67,77 @@ function test_nlexpr_scalarquadratic() return _test_nlexpr(f, terms, Dict(x => 1.1), 3.0) end +function test_nlexpr_scalarquadratic_1b() + x = MOI.VariableIndex(1) + f = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(1.1, x)], + [MOI.ScalarQuadraticTerm(2.5, x, x)], + 3.0, + ) + terms = [NL.OPMULT, 1.25, NL.OPMULT, x, x] + return _test_nlexpr(f, terms, Dict(x => 1.1), 3.0) +end + +function test_nlexpr_scalarquadratic_1c() + x = MOI.VariableIndex(1) + y = MOI.VariableIndex(2) + f = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(1.1, x)], + [MOI.ScalarQuadraticTerm(1.0, x, y)], + 3.0, + ) + terms = [NL.OPMULT, x, y] + return _test_nlexpr(f, terms, Dict(x => 1.1, y => 0.0), 3.0) +end + +function test_nlexpr_scalarquadratic_2() + x = MOI.VariableIndex(1) + y = MOI.VariableIndex(2) + f = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(1.1, x)], + [ + MOI.ScalarQuadraticTerm(2.0, x, x), + MOI.ScalarQuadraticTerm(1.0, x, y), + ], + 3.0, + ) + terms = [NL.OPPLUS, NL.OPMULT, x, x, NL.OPMULT, x, y] + return _test_nlexpr(f, terms, Dict(x => 1.1, y => 0.0), 3.0) +end + +function test_nlexpr_scalarquadratic_3() + x = MOI.VariableIndex(1) + y = MOI.VariableIndex(2) + z = MOI.VariableIndex(3) + f = MOI.ScalarQuadraticFunction( + [MOI.ScalarAffineTerm(1.1, x)], + [ + MOI.ScalarQuadraticTerm(2.0, x, x), + MOI.ScalarQuadraticTerm(0.5, x, y), + MOI.ScalarQuadraticTerm(4.0, x, z), + ], + 3.0, + ) + terms = [ + NL.OPSUMLIST, + 3, + NL.OPMULT, + x, + x, + NL.OPMULT, + 0.5, + NL.OPMULT, + x, + y, + NL.OPMULT, + 4.0, + NL.OPMULT, + x, + z, + ] + return _test_nlexpr(f, terms, Dict(x => 1.1, y => 0.0, z => 0.0), 3.0) +end + function test_nlexpr_unary_addition() x = MOI.VariableIndex(1) return _test_nlexpr(:(+$x), [x], Dict(x => 0), 0.0) @@ -121,6 +221,16 @@ function test_nlexpr_ref() return _test_nlexpr(:(x[$x]), [x], Dict(x => 0.0), 0.0) end +function test_nlexpr_empty() + return _test_nlexpr( + :(), + NL._NLTerm[], + Dict{MOI.VariableIndex,Float64}(), + 0.0; + force_nonlinear = true, + ) +end + function test_nlconstraint_interval() x = MOI.VariableIndex(1) expr = :(1.0 <= $x <= 2.0) @@ -874,9 +984,9 @@ function test_empty() y = MOI.add_variables(model, 3) n = NL.Model() map = MOI.copy_to(n, model) - @test !MOI.is_empty(model) - MOI.empty!(model) - @test MOI.is_empty(model) + @test !MOI.is_empty(n) + MOI.empty!(n) + @test MOI.is_empty(n) end function runtests() From 432a4e108138804c2d59ce3730476c0dc7ab4017 Mon Sep 17 00:00:00 2001 From: odow Date: Mon, 10 May 2021 12:49:26 +1200 Subject: [PATCH 10/12] Get to 100% coverage --- src/FileFormats/NL/NL.jl | 31 +++++++------------------------ src/FileFormats/NL/gen.jl | 20 -------------------- src/FileFormats/NL/opcode.jl | 13 ++++++++++++- test/FileFormats/NL/NL.jl | 13 +++++++++++++ 4 files changed, 32 insertions(+), 45 deletions(-) delete mode 100644 src/FileFormats/NL/gen.jl diff --git a/src/FileFormats/NL/NL.jl b/src/FileFormats/NL/NL.jl index a5fe132518..302b5a8f7f 100644 --- a/src/FileFormats/NL/NL.jl +++ b/src/FileFormats/NL/NL.jl @@ -445,14 +445,15 @@ function _process_constraint( op, l, u = _set_to_bounds(s) con = _NLConstraint(l, u, op, _NLExpr(f)) if isempty(con.expr.linear_terms) && isempty(con.expr.nonlinear_terms) - if l <= con.expr.constant <= u - continue - else + if !(l <= con.expr.constant <= u) error( "Malformed constraint. There are no variables and the " * "function constant $(con.expr.constant) is not in [$l, $u]", ) end + # Just use a placeholder for the constraint index. It's not going to + # be used. + mapping[ci] = MOI.ConstraintIndex{F,S}(-abs(ci.value)) elseif con.expr.is_linear push!(dest.h, con) mapping[ci] = MOI.ConstraintIndex{F,S}(length(dest.h)) @@ -514,31 +515,13 @@ function _write_nlexpr(io::IO, expr::_NLExpr, nlmodel::Model) _write_term(io, expr.constant, nlmodel) return end - # If the constant term is non-zero, we need to write it out. - skip_terms = 0 if !iszero(expr.constant) - if expr.nonlinear_terms[1] == OPSUMLIST - # The nonlinear expression is a summation. We can write our constant - # first, but we also need to increment the number of arguments by - # one. In addition, since we're writing out the first two terms now, - # we must skip them below. - _write_term(io, OPSUMLIST, nlmodel) - println(io, expr.nonlinear_terms[2] + 1) - _write_term(io, expr.constant, nlmodel) - skip_terms = 2 - else - # The nonlinear expression is something other than a summation, so - # add a new + node to the expression. - _write_term(io, OPPLUS, nlmodel) - _write_term(io, expr.constant, nlmodel) - end + # If the constant term is non-zero, we need to write it out. + _write_term(io, OPPLUS, nlmodel) + _write_term(io, expr.constant, nlmodel) end last_nary = false for term in expr.nonlinear_terms - if skip_terms > 0 - skip_terms -= 1 - continue - end if last_nary println(io, term::Int) last_nary = false diff --git a/src/FileFormats/NL/gen.jl b/src/FileFormats/NL/gen.jl deleted file mode 100644 index ed3321b864..0000000000 --- a/src/FileFormats/NL/gen.jl +++ /dev/null @@ -1,20 +0,0 @@ -# This script builds the list of recognized ASL opcodes using the header files -# in ASL_jll. -# -# Because it is used for code generation, it is not included by -# `MathOptInterface.FileFormats.NL`. Only re-run it if AMPL adds new opcodes -# (which is unlikely). - -using ASL_jll - -open("opcode.jl", "w") do io - println( - "# Do not modify. This file is automatically created by the script " * - "in `gen.jl`.\n", - ) - filename = joinpath(ASL_jll.artifact_dir, "include", "opcode.hd") - for line in readlines(filename) - items = split(line, c -> c == '\t' || c == ' '; keepempty = false) - println(io, "const ", items[2], " = ", items[3]) - end -end diff --git a/src/FileFormats/NL/opcode.jl b/src/FileFormats/NL/opcode.jl index 2a066b14ab..d62b9a3427 100644 --- a/src/FileFormats/NL/opcode.jl +++ b/src/FileFormats/NL/opcode.jl @@ -1,4 +1,15 @@ -# Do not modify. This file is automatically created by the script in `gen.jl`. +# This file is automatically created by the folowing script. Only re-run it if +# AMPL adds new opcodes (which is unlikely). +# ```julia +# using ASL_jll +# open("opcode.jl", "w") do io +# filename = joinpath(ASL_jll.artifact_dir, "include", "opcode.hd") +# for line in readlines(filename) +# items = split(line, c -> c == '\t' || c == ' '; keepempty = false) +# println(io, "const ", items[2], " = ", items[3]) +# end +# end +# ``` const OPPLUS = 0 const OPMINUS = 1 diff --git a/test/FileFormats/NL/NL.jl b/test/FileFormats/NL/NL.jl index ef0decb596..df0184c9bb 100644 --- a/test/FileFormats/NL/NL.jl +++ b/test/FileFormats/NL/NL.jl @@ -837,6 +837,19 @@ function test_issue_79() @test n.x[z].order == 0 end +function test_malformed_constraint_error() + model = MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()) + x = MOI.add_variable(model) + ci = MOI.add_constraint( + model, + MOI.ScalarAffineFunction(MOI.ScalarAffineTerm{Float64}[], 1.0), + MOI.LessThan(2.0), + ) + n = NL.Model() + mapping = MOI.copy_to(n, model) + @test haskey(mapping, ci) +end + function test_malformed_constraint_error() model = MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()) x = MOI.add_variable(model) From 8fafd20bc6bc673dbc1e06170b414399c6985310 Mon Sep 17 00:00:00 2001 From: odow Date: Mon, 10 May 2021 14:40:26 +1200 Subject: [PATCH 11/12] Fix test --- test/FileFormats/NL/NL.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/FileFormats/NL/NL.jl b/test/FileFormats/NL/NL.jl index df0184c9bb..69e05fca32 100644 --- a/test/FileFormats/NL/NL.jl +++ b/test/FileFormats/NL/NL.jl @@ -837,7 +837,7 @@ function test_issue_79() @test n.x[z].order == 0 end -function test_malformed_constraint_error() +function test_empty_constraint() model = MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()) x = MOI.add_variable(model) ci = MOI.add_constraint( From 454340d0767b6a8b4c893c3c3465f7978353cd33 Mon Sep 17 00:00:00 2001 From: odow Date: Mon, 10 May 2021 15:47:40 +1200 Subject: [PATCH 12/12] Revert mapping.conmap[F,S] change --- src/FileFormats/NL/NL.jl | 18 ++++++++++++------ 1 file changed, 12 insertions(+), 6 deletions(-) diff --git a/src/FileFormats/NL/NL.jl b/src/FileFormats/NL/NL.jl index 302b5a8f7f..66b96469c2 100644 --- a/src/FileFormats/NL/NL.jl +++ b/src/FileFormats/NL/NL.jl @@ -288,7 +288,7 @@ function MOI.copy_to( resize!(dest.order, length(dest.x)) # Now deal with the normal MOI constraints. for (F, S) in MOI.get(model, MOI.ListOfConstraintTypesPresent()) - _process_constraint(dest, model, mapping.conmap[F, S]) + _process_constraint(dest, model, F, S, mapping) end # Correct bounds of binary variables. Mainly because AMPL doesn't have the # concept of binary nonlinear variables, but it does have binary linear @@ -437,7 +437,9 @@ _set_to_bounds(set::MOI.EqualTo) = (4, set.value, set.value) function _process_constraint( dest::Model, model, - mapping::MOI.Utilities.DoubleDicts.IndexWithType{F,S}, + ::Type{F}, + ::Type{S}, + mapping, ) where {F,S} for ci in MOI.get(model, MOI.ListOfConstraintIndices{F,S}()) f = MOI.get(model, MOI.ConstraintFunction(), ci) @@ -468,8 +470,10 @@ end function _process_constraint( dest::Model, model, - mapping::MOI.Utilities.DoubleDicts.IndexWithType{F,S}, -) where {F<:MOI.SingleVariable,S} + F::Type{MOI.SingleVariable}, + S::Type{<:_SCALAR_SETS}, + mapping, +) for ci in MOI.get(model, MOI.ListOfConstraintIndices{F,S}()) mapping[ci] = ci f = MOI.get(model, MOI.ConstraintFunction(), ci) @@ -488,8 +492,10 @@ end function _process_constraint( dest::Model, model, - mapping::MOI.Utilities.DoubleDicts.IndexWithType{F,S}, -) where {F<:MOI.SingleVariable,S<:Union{MOI.ZeroOne,MOI.Integer}} + F::Type{MOI.SingleVariable}, + S::Type{<:Union{MOI.ZeroOne,MOI.Integer}}, + mapping, +) for ci in MOI.get(model, MOI.ListOfConstraintIndices{F,S}()) mapping[ci] = ci f = MOI.get(model, MOI.ConstraintFunction(), ci)