From 4b44d124c00d52e99ec2b0f8605d8a5c285c117c Mon Sep 17 00:00:00 2001 From: Vaibhav Dixit Date: Thu, 23 Nov 2023 00:31:10 -0500 Subject: [PATCH 1/7] Keep symbolic expressions from mtk, make it moi compatible later --- lib/OptimizationMOI/src/OptimizationMOI.jl | 21 +++++++++++++++++ lib/OptimizationMOI/src/moi.jl | 27 ++++------------------ lib/OptimizationMOI/src/nlp.jl | 18 ++++++++++++--- 3 files changed, 40 insertions(+), 26 deletions(-) diff --git a/lib/OptimizationMOI/src/OptimizationMOI.jl b/lib/OptimizationMOI/src/OptimizationMOI.jl index 4d871b108..d8c814393 100644 --- a/lib/OptimizationMOI/src/OptimizationMOI.jl +++ b/lib/OptimizationMOI/src/OptimizationMOI.jl @@ -116,6 +116,27 @@ function __moi_status_to_ReturnCode(status::MOI.TerminationStatusCode) end end +function get_expr_map(prob, f) + pairs_arr = if prob.p isa SciMLBase.NullParameters + [_s => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)] + else + vcat([_s => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)], + [_p => p[i] for (i, _p) in enumerate(f.paramsyms)]) + end + + return pairs_arr +end + +""" +Substitute variables and parameters with canonical names x and p respectively. +""" +function rep_pars_vals!(e::Expr, p) + rep_pars_vals!.(e.args, Ref(p)) + replace!(e.args, p...) +end + +function rep_pars_vals!(e, p) end + """ Replaces every expression `:x[i]` with `:x[MOI.VariableIndex(i)]` """ diff --git a/lib/OptimizationMOI/src/moi.jl b/lib/OptimizationMOI/src/moi.jl index 3546637a4..accdce4e8 100644 --- a/lib/OptimizationMOI/src/moi.jl +++ b/lib/OptimizationMOI/src/moi.jl @@ -17,17 +17,13 @@ function MOIOptimizationCache(prob::OptimizationProblem, opt; kwargs...) throw(ArgumentError("Expected an `OptimizationProblem` that was setup via an `OptimizationSystem`, consider `modelingtoolkitize(prob).`")) # TODO: check if the problem is at most bilinear, i.e. affine and or quadratic terms in two variables - expr_map = get_expr_map(prob.f.sys) - expr = repl_getindex!(convert_to_expr(MTK.subs_constants(MTK.objective(prob.f.sys)), - prob.f.sys; expand_expr = false, expr_map)) + expr_map = get_expr_map(prob, prob.f) + expr = repl_getindex!(convert_to_expr(prob.f.obj_expr, expr_map; expand_expr = false)) cons = MTK.constraints(prob.f.sys) cons_expr = Vector{Expr}(undef, length(cons)) Threads.@sync for i in eachindex(cons) - Threads.@spawn cons_expr[i] = repl_getindex!(convert_to_expr(Symbolics.canonical_form(MTK.subs_constants(cons[i])), - prob.f.sys; - expand_expr = false, - expr_map)) + Threads.@spawn cons_expr[i] = repl_getindex!(convert_to_expr(prob.f.cons_expr, expr_map; expand_expr = false)) end return MOIOptimizationCache(prob.f, @@ -421,21 +417,6 @@ function symbolify!(e::Expr) return e end -""" - get_expr_map(sys) - -Make a map from every parameter and state of the given system to an expression indexing its position -in the state or parameter vector. -""" -function get_expr_map(sys) - dvs = ModelingToolkit.states(sys) - ps = ModelingToolkit.parameters(sys) - return vcat([ModelingToolkit.toexpr(_s) => Expr(:ref, :x, i) - for (i, _s) in enumerate(dvs)], - [ModelingToolkit.toexpr(_p) => Expr(:ref, :p, i) - for (i, _p) in enumerate(ps)]) -end - """ convert_to_expr(eq, sys; expand_expr = false, pairs_arr = expr_map(sys)) @@ -447,7 +428,7 @@ parameters with `x[i]` and `p[i]`. - `sys`: Reference to the system holding the parameters and states - `expand_expr=false`: If `true` the symbolic expression is expanded first. """ -function convert_to_expr(eq, sys; expand_expr = false, expr_map = get_expr_map(sys)) +function convert_to_expr(eq, expr_map; expand_expr = false) if expand_expr eq = try Symbolics.expand(eq) # PolyForm sometimes errors diff --git a/lib/OptimizationMOI/src/nlp.jl b/lib/OptimizationMOI/src/nlp.jl index e1c30d40e..ad1635f53 100644 --- a/lib/OptimizationMOI/src/nlp.jl +++ b/lib/OptimizationMOI/src/nlp.jl @@ -15,6 +15,8 @@ mutable struct MOIOptimizationNLPEvaluator{T, F <: OptimizationFunction, RC, LB, H::HT cons_H::Vector{CHT} callback::CB + obj_expr::Union{Expr, Nothing} + cons_expr::Union{Vector{Expr}, Nothing} end function Base.getproperty(evaluator::MOIOptimizationNLPEvaluator, x::Symbol) @@ -135,6 +137,14 @@ function MOIOptimizationNLPCache(prob::OptimizationProblem, end lcons = prob.lcons === nothing ? fill(T(-Inf), num_cons) : prob.lcons ucons = prob.ucons === nothing ? fill(T(Inf), num_cons) : prob.ucons + + obj_expr = toexpr(f.expr) + cons_expr = toexpr.(f.cons_expr) + + pairs_arr = get_expr_map(prob, f) + + rep_pars_vals!(obj_expr, pairs_arr) + rep_pars_vals!.(cons_expr, Ref(pairs_arr)) evaluator = MOIOptimizationNLPEvaluator(f, reinit_cache, @@ -147,7 +157,9 @@ function MOIOptimizationNLPCache(prob::OptimizationProblem, J, H, cons_H, - callback) + callback, + obj_expr, + cons_expr) return MOIOptimizationNLPCache(evaluator, opt, NamedTuple(kwargs)) end @@ -334,7 +346,7 @@ function MOI.eval_hessian_lagrangian(evaluator::MOIOptimizationNLPEvaluator{T}, end function MOI.objective_expr(evaluator::MOIOptimizationNLPEvaluator) - expr = deepcopy(evaluator.f.expr) + expr = deepcopy(evaluator.obj_expr) repl_getindex!(expr) _replace_parameter_indices!(expr, evaluator.p) _replace_variable_indices!(expr) @@ -343,7 +355,7 @@ end function MOI.constraint_expr(evaluator::MOIOptimizationNLPEvaluator, i) # expr has the form f(x,p) == 0 or f(x,p) <= 0 - cons_expr = deepcopy(evaluator.f.cons_expr[i].args[2]) + cons_expr = deepcopy(evaluator.cons_expr[i].args[2]) repl_getindex!(cons_expr) _replace_parameter_indices!(cons_expr, evaluator.p) _replace_variable_indices!(cons_expr) From 822b67784415664afb35a28e24a2923c98b91f22 Mon Sep 17 00:00:00 2001 From: Vaibhav Dixit Date: Mon, 27 Nov 2023 12:55:34 -0500 Subject: [PATCH 2/7] wip --- lib/OptimizationMOI/src/OptimizationMOI.jl | 7 +++--- lib/OptimizationMOI/src/moi.jl | 25 ++++++++++++++-------- 2 files changed, 19 insertions(+), 13 deletions(-) diff --git a/lib/OptimizationMOI/src/OptimizationMOI.jl b/lib/OptimizationMOI/src/OptimizationMOI.jl index d8c814393..06b48882f 100644 --- a/lib/OptimizationMOI/src/OptimizationMOI.jl +++ b/lib/OptimizationMOI/src/OptimizationMOI.jl @@ -118,12 +118,11 @@ end function get_expr_map(prob, f) pairs_arr = if prob.p isa SciMLBase.NullParameters - [_s => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)] + [Meta.parse(string(_s)) => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)] else - vcat([_s => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)], - [_p => p[i] for (i, _p) in enumerate(f.paramsyms)]) + vcat([Symbol(_s) => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)], + [Symbol(_s) => Expr(:ref, p, i) for (i, _p) in enumerate(f.paramsyms)]) end - return pairs_arr end diff --git a/lib/OptimizationMOI/src/moi.jl b/lib/OptimizationMOI/src/moi.jl index accdce4e8..457ab49be 100644 --- a/lib/OptimizationMOI/src/moi.jl +++ b/lib/OptimizationMOI/src/moi.jl @@ -13,21 +13,28 @@ struct MOIOptimizationCache{F <: OptimizationFunction, RC, LB, UB, I, S, EX, end function MOIOptimizationCache(prob::OptimizationProblem, opt; kwargs...) - isnothing(prob.f.sys) && - throw(ArgumentError("Expected an `OptimizationProblem` that was setup via an `OptimizationSystem`, consider `modelingtoolkitize(prob).`")) - + f = prob.f + reinit_cache = Optimization.ReInitCache(prob.u0, prob.p) + if isnothing(f.sys) + if f.adtype isa Optimization.AutoModelingToolkit + num_cons = prob.ucons === nothing ? 0 : length(prob.ucons) + f = Optimization.instatiate_function(prob.f, reinit_cache, prob.f.adtype, num_cons) + else + throw(ArgumentError("Expected an `OptimizationProblem` that was setup via an `OptimizationSystem`, or AutoModelingToolkit ad choice")) + end + end # TODO: check if the problem is at most bilinear, i.e. affine and or quadratic terms in two variables - expr_map = get_expr_map(prob, prob.f) - expr = repl_getindex!(convert_to_expr(prob.f.obj_expr, expr_map; expand_expr = false)) + expr_map = get_expr_map(prob, f) + expr = repl_getindex!(convert_to_expr(f.expr, expr_map; expand_expr = false)) - cons = MTK.constraints(prob.f.sys) + cons = MTK.constraints(f.sys) cons_expr = Vector{Expr}(undef, length(cons)) Threads.@sync for i in eachindex(cons) - Threads.@spawn cons_expr[i] = repl_getindex!(convert_to_expr(prob.f.cons_expr, expr_map; expand_expr = false)) + Threads.@spawn cons_expr[i] = repl_getindex!(convert_to_expr(f.cons_expr[i], expr_map; expand_expr = false)) end - return MOIOptimizationCache(prob.f, - Optimization.ReInitCache(prob.u0, prob.p), + return MOIOptimizationCache(f, + reinit_cache, prob.lb, prob.ub, prob.int, From 8ab2e6fa109a7969ec774032a5d10cf48111e820 Mon Sep 17 00:00:00 2001 From: Vaibhav Dixit Date: Sat, 2 Dec 2023 11:36:03 -0500 Subject: [PATCH 3/7] wip --- lib/OptimizationMOI/src/OptimizationMOI.jl | 34 +++++++++++----------- lib/OptimizationMOI/src/moi.jl | 12 +++++++- 2 files changed, 28 insertions(+), 18 deletions(-) diff --git a/lib/OptimizationMOI/src/OptimizationMOI.jl b/lib/OptimizationMOI/src/OptimizationMOI.jl index 06b48882f..1f3d3a66d 100644 --- a/lib/OptimizationMOI/src/OptimizationMOI.jl +++ b/lib/OptimizationMOI/src/OptimizationMOI.jl @@ -116,25 +116,25 @@ function __moi_status_to_ReturnCode(status::MOI.TerminationStatusCode) end end -function get_expr_map(prob, f) - pairs_arr = if prob.p isa SciMLBase.NullParameters - [Meta.parse(string(_s)) => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)] - else - vcat([Symbol(_s) => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)], - [Symbol(_s) => Expr(:ref, p, i) for (i, _p) in enumerate(f.paramsyms)]) - end - return pairs_arr -end +# function get_expr_map(prob, f) +# pairs_arr = if prob.p isa SciMLBase.NullParameters +# [Meta.parse(string(_s)) => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)] +# else +# vcat([Symbol(_s) => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)], +# [Symbol(_s) => Expr(:ref, p, i) for (i, _p) in enumerate(f.paramsyms)]) +# end +# return pairs_arr +# end -""" -Substitute variables and parameters with canonical names x and p respectively. -""" -function rep_pars_vals!(e::Expr, p) - rep_pars_vals!.(e.args, Ref(p)) - replace!(e.args, p...) -end +# """ +# Substitute variables and parameters with canonical names x and p respectively. +# """ +# function rep_pars_vals!(e::Expr, p) +# rep_pars_vals!.(e.args, Ref(p)) +# replace!(e.args, p...) +# end -function rep_pars_vals!(e, p) end +# function rep_pars_vals!(e, p) end """ Replaces every expression `:x[i]` with `:x[MOI.VariableIndex(i)]` diff --git a/lib/OptimizationMOI/src/moi.jl b/lib/OptimizationMOI/src/moi.jl index 457ab49be..f6414ac54 100644 --- a/lib/OptimizationMOI/src/moi.jl +++ b/lib/OptimizationMOI/src/moi.jl @@ -23,8 +23,9 @@ function MOIOptimizationCache(prob::OptimizationProblem, opt; kwargs...) throw(ArgumentError("Expected an `OptimizationProblem` that was setup via an `OptimizationSystem`, or AutoModelingToolkit ad choice")) end end + # TODO: check if the problem is at most bilinear, i.e. affine and or quadratic terms in two variables - expr_map = get_expr_map(prob, f) + expr_map = get_expr_map(prob.f.sys) expr = repl_getindex!(convert_to_expr(f.expr, expr_map; expand_expr = false)) cons = MTK.constraints(f.sys) @@ -448,3 +449,12 @@ function convert_to_expr(eq, expr_map; expand_expr = false) expr = symbolify!(expr) return expr end + +function get_expr_map(sys) + dvs = ModelingToolkit.states(sys) + ps = ModelingToolkit.parameters(sys) + return vcat([ModelingToolkit.toexpr(_s) => Expr(:ref, :x, i) + for (i, _s) in enumerate(dvs)], + [ModelingToolkit.toexpr(_p) => Expr(:ref, :p, i) + for (i, _p) in enumerate(ps)]) +end \ No newline at end of file From 9f568c26bb848cd7a44328e58906f6bdd04e975d Mon Sep 17 00:00:00 2001 From: Vaibhav Dixit Date: Sat, 16 Dec 2023 23:01:43 -0500 Subject: [PATCH 4/7] using Symbolics --- lib/OptimizationMOI/src/OptimizationMOI.jl | 2 +- lib/OptimizationMOI/src/moi.jl | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/lib/OptimizationMOI/src/OptimizationMOI.jl b/lib/OptimizationMOI/src/OptimizationMOI.jl index 1f3d3a66d..cb2cce88d 100644 --- a/lib/OptimizationMOI/src/OptimizationMOI.jl +++ b/lib/OptimizationMOI/src/OptimizationMOI.jl @@ -8,7 +8,7 @@ using SparseArrays import ModelingToolkit using ModelingToolkit: parameters, states, varmap_to_vars, mergedefaults const MTK = ModelingToolkit -import Symbolics +using Symbolics const MOI = MathOptInterface diff --git a/lib/OptimizationMOI/src/moi.jl b/lib/OptimizationMOI/src/moi.jl index f6414ac54..63cc9e62c 100644 --- a/lib/OptimizationMOI/src/moi.jl +++ b/lib/OptimizationMOI/src/moi.jl @@ -18,7 +18,7 @@ function MOIOptimizationCache(prob::OptimizationProblem, opt; kwargs...) if isnothing(f.sys) if f.adtype isa Optimization.AutoModelingToolkit num_cons = prob.ucons === nothing ? 0 : length(prob.ucons) - f = Optimization.instatiate_function(prob.f, reinit_cache, prob.f.adtype, num_cons) + f = Optimization.instantiate_function(prob.f, reinit_cache, prob.f.adtype, num_cons) else throw(ArgumentError("Expected an `OptimizationProblem` that was setup via an `OptimizationSystem`, or AutoModelingToolkit ad choice")) end From 6a1f353007df97d16853679426c7a1ea2a011a77 Mon Sep 17 00:00:00 2001 From: Vaibhav Dixit Date: Sat, 16 Dec 2023 23:28:15 -0500 Subject: [PATCH 5/7] MTK toexpr? --- lib/OptimizationMOI/src/OptimizationMOI.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/lib/OptimizationMOI/src/OptimizationMOI.jl b/lib/OptimizationMOI/src/OptimizationMOI.jl index cb2cce88d..863998c93 100644 --- a/lib/OptimizationMOI/src/OptimizationMOI.jl +++ b/lib/OptimizationMOI/src/OptimizationMOI.jl @@ -5,8 +5,8 @@ using Reexport using MathOptInterface using Optimization.SciMLBase using SparseArrays +import ModelingToolkit: parameters, states, varmap_to_vars, mergedefaults, toexpr import ModelingToolkit -using ModelingToolkit: parameters, states, varmap_to_vars, mergedefaults const MTK = ModelingToolkit using Symbolics From 08af6c00537e4fa945ecf64804c4d8c358eff3de Mon Sep 17 00:00:00 2001 From: Vaibhav Dixit Date: Sun, 17 Dec 2023 23:26:15 -0500 Subject: [PATCH 6/7] Everything other than amplnlwriter tests working --- lib/OptimizationMOI/Project.toml | 2 + lib/OptimizationMOI/src/OptimizationMOI.jl | 116 ++++++++++++++++++--- lib/OptimizationMOI/src/moi.jl | 106 +------------------ lib/OptimizationMOI/src/nlp.jl | 39 +++++-- 4 files changed, 134 insertions(+), 129 deletions(-) diff --git a/lib/OptimizationMOI/Project.toml b/lib/OptimizationMOI/Project.toml index aa504e02a..a5666dc7a 100644 --- a/lib/OptimizationMOI/Project.toml +++ b/lib/OptimizationMOI/Project.toml @@ -4,12 +4,14 @@ authors = ["Vaibhav Dixit and contributors"] version = "0.1.16" [deps] +AmplNLWriter = "7c4d4715-977e-5154-bfe0-e096adeac482" Ipopt_jll = "9cc047cb-c261-5740-88fc-0cf96f7bdcc7" MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78" Optimization = "7f7a1694-90dd-40f0-9382-eb1efda571ba" Reexport = "189a3867-3050-52da-a836-e630ba90ab69" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" +SymbolicIndexingInterface = "2efcf032-c050-4f8e-a9bb-153293bab1f5" Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7" [compat] diff --git a/lib/OptimizationMOI/src/OptimizationMOI.jl b/lib/OptimizationMOI/src/OptimizationMOI.jl index 863998c93..fd3b8176f 100644 --- a/lib/OptimizationMOI/src/OptimizationMOI.jl +++ b/lib/OptimizationMOI/src/OptimizationMOI.jl @@ -4,6 +4,7 @@ using Reexport @reexport using Optimization using MathOptInterface using Optimization.SciMLBase +using SymbolicIndexingInterface using SparseArrays import ModelingToolkit: parameters, states, varmap_to_vars, mergedefaults, toexpr import ModelingToolkit @@ -116,25 +117,106 @@ function __moi_status_to_ReturnCode(status::MOI.TerminationStatusCode) end end -# function get_expr_map(prob, f) -# pairs_arr = if prob.p isa SciMLBase.NullParameters -# [Meta.parse(string(_s)) => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)] -# else -# vcat([Symbol(_s) => Expr(:ref, :x, i) for (i, _s) in enumerate(f.syms)], -# [Symbol(_s) => Expr(:ref, p, i) for (i, _p) in enumerate(f.paramsyms)]) -# end -# return pairs_arr -# end +_get_variable_index_from_expr(expr::T) where {T} = throw(MalformedExprException("$expr")) +function _get_variable_index_from_expr(expr::Expr) + _is_var_ref_expr(expr) + return MOI.VariableIndex(expr.args[2]) +end -# """ -# Substitute variables and parameters with canonical names x and p respectively. -# """ -# function rep_pars_vals!(e::Expr, p) -# rep_pars_vals!.(e.args, Ref(p)) -# replace!(e.args, p...) -# end +function _is_var_ref_expr(expr::Expr) + expr.head == :ref || throw(MalformedExprException("$expr")) # x[i] + expr.args[1] == :x || throw(MalformedExprException("$expr")) + return true +end + +function is_eq(expr::Expr) + expr.head == :call || throw(MalformedExprException("$expr")) + expr.args[1] in [:(==), :(=)] +end + +function is_leq(expr::Expr) + expr.head == :call || throw(MalformedExprException("$expr")) + expr.args[1] == :(<=) +end + +""" + rep_pars_vals!(expr::T, expr_map) + +Replaces variable expressions of the form `:some_variable` or `:(getindex, :some_variable, j)` with +`x[i]` were `i` is the corresponding index in the state vector. Same for the parameters. The +variable/parameter pairs are provided via the `expr_map`. + +Expects only expressions where the variables and parameters are of the form `:some_variable` +or `:(getindex, :some_variable, j)` or :(some_variable[j]). +""" +rep_pars_vals!(expr::T, expr_map) where {T} = expr +function rep_pars_vals!(expr::Symbol, expr_map) + for (f, n) in expr_map + isequal(f, expr) && return n + end + return expr +end +function rep_pars_vals!(expr::Expr, expr_map) + if (expr.head == :call && expr.args[1] == getindex) || (expr.head == :ref) + for (f, n) in expr_map + isequal(f, expr) && return n + end + end + Threads.@sync for i in eachindex(expr.args) + i == 1 && expr.head == :call && continue # first arg is the operator + Threads.@spawn expr.args[i] = rep_pars_vals!(expr.args[i], expr_map) + end + return expr +end + +""" + symbolify!(e) + +Ensures that a given expression is fully symbolic, e.g. no function calls. +""" +symbolify!(e) = e +function symbolify!(e::Expr) + if !(e.args[1] isa Symbol) + e.args[1] = Symbol(e.args[1]) + end + symbolify!.(e.args) + return e +end + +""" + convert_to_expr(eq, sys; expand_expr = false, pairs_arr = expr_map(sys)) + +Converts the given symbolic expression to a Julia `Expr` and replaces all symbols, i.e. states and +parameters with `x[i]` and `p[i]`. + +# Arguments: +- `eq`: Expression to convert +- `sys`: Reference to the system holding the parameters and states +- `expand_expr=false`: If `true` the symbolic expression is expanded first. +""" +function convert_to_expr(eq, expr_map; expand_expr = false) + if expand_expr + eq = try + Symbolics.expand(eq) # PolyForm sometimes errors + catch e + Symbolics.expand(eq) + end + end + expr = ModelingToolkit.toexpr(eq) -# function rep_pars_vals!(e, p) end + expr = rep_pars_vals!(expr, expr_map) + expr = symbolify!(expr) + return expr +end + +function get_expr_map(sys) + dvs = ModelingToolkit.states(sys) + ps = ModelingToolkit.parameters(sys) + return vcat([ModelingToolkit.toexpr(_s) => Expr(:ref, :x, i) + for (i, _s) in enumerate(dvs)], + [ModelingToolkit.toexpr(_p) => Expr(:ref, :p, i) + for (i, _p) in enumerate(ps)]) +end """ Replaces every expression `:x[i]` with `:x[MOI.VariableIndex(i)]` diff --git a/lib/OptimizationMOI/src/moi.jl b/lib/OptimizationMOI/src/moi.jl index 63cc9e62c..dba0a5207 100644 --- a/lib/OptimizationMOI/src/moi.jl +++ b/lib/OptimizationMOI/src/moi.jl @@ -26,8 +26,8 @@ function MOIOptimizationCache(prob::OptimizationProblem, opt; kwargs...) # TODO: check if the problem is at most bilinear, i.e. affine and or quadratic terms in two variables expr_map = get_expr_map(prob.f.sys) - expr = repl_getindex!(convert_to_expr(f.expr, expr_map; expand_expr = false)) - + expr = convert_to_expr(f.expr, expr_map; expand_expr = false) + expr = repl_getindex!(expr) cons = MTK.constraints(f.sys) cons_expr = Vector{Expr}(undef, length(cons)) Threads.@sync for i in eachindex(cons) @@ -356,105 +356,3 @@ function collect_moi_terms!(expr::Expr, affine_terms, quadratic_terms, constant) return end - -_get_variable_index_from_expr(expr::T) where {T} = throw(MalformedExprException("$expr")) -function _get_variable_index_from_expr(expr::Expr) - _is_var_ref_expr(expr) - return MOI.VariableIndex(expr.args[2]) -end - -function _is_var_ref_expr(expr::Expr) - expr.head == :ref || throw(MalformedExprException("$expr")) # x[i] - expr.args[1] == :x || throw(MalformedExprException("$expr")) - return true -end - -function is_eq(expr::Expr) - expr.head == :call || throw(MalformedExprException("$expr")) - expr.args[1] in [:(==), :(=)] -end - -function is_leq(expr::Expr) - expr.head == :call || throw(MalformedExprException("$expr")) - expr.args[1] == :(<=) -end - -############################################################################################## -## TODO: remove if in ModelingToolkit -""" - rep_pars_vals!(expr::T, expr_map) - -Replaces variable expressions of the form `:some_variable` or `:(getindex, :some_variable, j)` with -`x[i]` were `i` is the corresponding index in the state vector. Same for the parameters. The -variable/parameter pairs are provided via the `expr_map`. - -Expects only expressions where the variables and parameters are of the form `:some_variable` -or `:(getindex, :some_variable, j)` or :(some_variable[j]). -""" -rep_pars_vals!(expr::T, expr_map) where {T} = expr -function rep_pars_vals!(expr::Symbol, expr_map) - for (f, n) in expr_map - isequal(f, expr) && return n - end - return expr -end -function rep_pars_vals!(expr::Expr, expr_map) - if (expr.head == :call && expr.args[1] == getindex) || (expr.head == :ref) - for (f, n) in expr_map - isequal(f, expr) && return n - end - end - Threads.@sync for i in eachindex(expr.args) - i == 1 && expr.head == :call && continue # first arg is the operator - Threads.@spawn expr.args[i] = rep_pars_vals!(expr.args[i], expr_map) - end - return expr -end - -""" - symbolify!(e) - -Ensures that a given expression is fully symbolic, e.g. no function calls. -""" -symbolify!(e) = e -function symbolify!(e::Expr) - if !(e.args[1] isa Symbol) - e.args[1] = Symbol(e.args[1]) - end - symbolify!.(e.args) - return e -end - -""" - convert_to_expr(eq, sys; expand_expr = false, pairs_arr = expr_map(sys)) - -Converts the given symbolic expression to a Julia `Expr` and replaces all symbols, i.e. states and -parameters with `x[i]` and `p[i]`. - -# Arguments: -- `eq`: Expression to convert -- `sys`: Reference to the system holding the parameters and states -- `expand_expr=false`: If `true` the symbolic expression is expanded first. -""" -function convert_to_expr(eq, expr_map; expand_expr = false) - if expand_expr - eq = try - Symbolics.expand(eq) # PolyForm sometimes errors - catch e - Symbolics.expand(eq) - end - end - expr = ModelingToolkit.toexpr(eq) - expr = rep_pars_vals!(expr, expr_map) - expr = symbolify!(expr) - return expr -end - -function get_expr_map(sys) - dvs = ModelingToolkit.states(sys) - ps = ModelingToolkit.parameters(sys) - return vcat([ModelingToolkit.toexpr(_s) => Expr(:ref, :x, i) - for (i, _s) in enumerate(dvs)], - [ModelingToolkit.toexpr(_p) => Expr(:ref, :p, i) - for (i, _p) in enumerate(ps)]) -end \ No newline at end of file diff --git a/lib/OptimizationMOI/src/nlp.jl b/lib/OptimizationMOI/src/nlp.jl index ad1635f53..9a1bf5cbe 100644 --- a/lib/OptimizationMOI/src/nlp.jl +++ b/lib/OptimizationMOI/src/nlp.jl @@ -137,14 +137,37 @@ function MOIOptimizationNLPCache(prob::OptimizationProblem, end lcons = prob.lcons === nothing ? fill(T(-Inf), num_cons) : prob.lcons ucons = prob.ucons === nothing ? fill(T(Inf), num_cons) : prob.ucons - - obj_expr = toexpr(f.expr) - cons_expr = toexpr.(f.cons_expr) - pairs_arr = get_expr_map(prob, f) + if f.sys isa SymbolicIndexingInterface.SymbolCache{Nothing, Nothing, Nothing} + sys = MTK.modelingtoolkitize(prob) + if !isnothing(prob.p) && !(prob.p isa SciMLBase.NullParameters) + unames = variable_symbols(sys) + pnames = parameter_symbols(sys) + us = [unames[i] => prob.u0[i] for i in 1:length(prob.u0)] + ps = [pnames[i] => prob.p[i] for i in 1:length(prob.p)] + sysprob = OptimizationProblem(sys, us, ps) + else + unames = variable_symbols(sys) + us = [unames[i] => prob.u0[i] for i in 1:length(prob.u0)] + sysprob = OptimizationProblem(sys, us) + end - rep_pars_vals!(obj_expr, pairs_arr) - rep_pars_vals!.(cons_expr, Ref(pairs_arr)) + obj_expr = sysprob.f.expr + cons_expr = sysprob.f.cons_expr + else + sys = f.sys + obj_expr = f.expr + cons_expr = f.cons_expr + end + + expr_map = get_expr_map(sys) + expr = convert_to_expr(obj_expr, expr_map; expand_expr = false) + expr = repl_getindex!(expr) + cons = MTK.constraints(sys) + _cons_expr = Vector{Expr}(undef, length(cons)) + for i in eachindex(cons) + _cons_expr[i] = repl_getindex!(convert_to_expr(cons_expr[i], expr_map; expand_expr = false)) + end evaluator = MOIOptimizationNLPEvaluator(f, reinit_cache, @@ -158,8 +181,8 @@ function MOIOptimizationNLPCache(prob::OptimizationProblem, H, cons_H, callback, - obj_expr, - cons_expr) + expr, + _cons_expr) return MOIOptimizationNLPCache(evaluator, opt, NamedTuple(kwargs)) end From 75d197b71273c09880d5ecc830ee416bf36c43ab Mon Sep 17 00:00:00 2001 From: Vaibhav Dixit Date: Tue, 19 Dec 2023 13:40:41 -0500 Subject: [PATCH 7/7] Handle nlp constraints expressions --- lib/OptimizationMOI/Project.toml | 1 - lib/OptimizationMOI/src/nlp.jl | 6 +++--- 2 files changed, 3 insertions(+), 4 deletions(-) diff --git a/lib/OptimizationMOI/Project.toml b/lib/OptimizationMOI/Project.toml index a5666dc7a..52fb07202 100644 --- a/lib/OptimizationMOI/Project.toml +++ b/lib/OptimizationMOI/Project.toml @@ -4,7 +4,6 @@ authors = ["Vaibhav Dixit and contributors"] version = "0.1.16" [deps] -AmplNLWriter = "7c4d4715-977e-5154-bfe0-e096adeac482" Ipopt_jll = "9cc047cb-c261-5740-88fc-0cf96f7bdcc7" MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78" diff --git a/lib/OptimizationMOI/src/nlp.jl b/lib/OptimizationMOI/src/nlp.jl index 9a1bf5cbe..6b2289d8d 100644 --- a/lib/OptimizationMOI/src/nlp.jl +++ b/lib/OptimizationMOI/src/nlp.jl @@ -377,13 +377,13 @@ function MOI.objective_expr(evaluator::MOIOptimizationNLPEvaluator) end function MOI.constraint_expr(evaluator::MOIOptimizationNLPEvaluator, i) - # expr has the form f(x,p) == 0 or f(x,p) <= 0 + # expr has the form f(x,p) == 0 or f(x,p) <= 0 cons_expr = deepcopy(evaluator.cons_expr[i].args[2]) + compop = Symbol(evaluator.cons_expr[i].args[1]) repl_getindex!(cons_expr) _replace_parameter_indices!(cons_expr, evaluator.p) _replace_variable_indices!(cons_expr) - lb, ub = Float64(evaluator.lcons[i]), Float64(evaluator.ucons[i]) - return :($lb <= $cons_expr <= $ub) + return Expr(:call, compop, cons_expr, 0.0) end function _add_moi_variables!(opt_setup, evaluator::MOIOptimizationNLPEvaluator)