From 65eafbeef08be1b71f7380f20795b2b4c1703889 Mon Sep 17 00:00:00 2001 From: odow Date: Fri, 21 Apr 2023 17:39:19 +1200 Subject: [PATCH] WIP: update to new JuMP nonlinear syntax * All instances are a straight swap, dropping the NL part of macro names. * Duplicated code blocks that were added to work-around NL macros have been removed * Tests have been updated to new solution values Questions: * Should the test values have changed? Macro parsing might give slightly different numerical values. Would that make a difference? * I'll comment in-line, but are the higher-order polynomial objectives correct? --- Project.toml | 19 +++-- src/core/objective.jl | 42 +++++----- src/form/acp.jl | 187 +++++++++++++++--------------------------- src/form/acr.jl | 68 +++++---------- src/form/act.jl | 2 +- src/form/iv.jl | 44 +++++----- test/docs.jl | 2 +- test/model.jl | 2 +- test/opf-var.jl | 8 +- test/opf.jl | 48 ++++++----- test/runtests.jl | 13 +-- 11 files changed, 178 insertions(+), 257 deletions(-) diff --git a/Project.toml b/Project.toml index 51fd54c60..d0fc875b1 100644 --- a/Project.toml +++ b/Project.toml @@ -11,18 +11,19 @@ JuMP = "4076af6c-e467-56ae-b986-b466b2749572" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Memento = "f28f55f0-a522-5efc-85c2-fe41dfb9b2d9" NLsolve = "2774e3e8-f4cf-5e23-947b-6d7e65073b56" +Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" [compat] -HiGHS = "~0.3, ~1" -InfrastructureModels = "~0.6, ~0.7" -Ipopt = "~0.8, ~0.9, ~1" -JSON = "~0.18, ~0.19, ~0.20, ~0.21" -JuMP = "~0.22, ~0.23, 1" -Juniper = "~0.8, ~0.9" -Memento = "~1.0, ~1.1, ~1.2, ~1.3, ~1.4" -NLsolve = "4.0" -SCS = "~0.9, ~1.0" +HiGHS = "1" +InfrastructureModels = "0.6, 0.7" +Ipopt = "1" +JSON = "0.18, 0.19, 0.20, 0.21" +JuMP = "1" +Juniper = "0.8, 0.9" +Memento = "1" +NLsolve = "4" +SCS = "1" julia = "1.6" [extras] diff --git a/src/core/objective.jl b/src/core/objective.jl index b44b121c4..acf245ae0 100644 --- a/src/core/objective.jl +++ b/src/core/objective.jl @@ -276,16 +276,17 @@ function _objective_min_fuel_and_flow_cost_polynomial_nl(pm::AbstractPowerModel; cost_rev = reverse(gen["cost"]) if length(cost_rev) == 1 - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1]) + gen_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1]) elseif length(cost_rev) == 2 - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1] + cost_rev[2]*pg) + gen_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1] + cost_rev[2]*pg) elseif length(cost_rev) == 3 - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1] + cost_rev[2]*pg + cost_rev[3]*pg^2) + gen_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1] + cost_rev[2]*pg + cost_rev[3]*pg^2) elseif length(cost_rev) >= 4 cost_rev_nl = cost_rev[4:end] - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1] + cost_rev[2]*pg + cost_rev[3]*pg^2 + sum( v*pg^(d+3) for (d,v) in enumerate(cost_rev_nl)) ) + # TODO(odow): is the correct? Or should it be pg^(d+2.0)? + gen_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1] + cost_rev[2]*pg + cost_rev[3]*pg^2 + sum( v*pg^(d+3.0) for (d,v) in enumerate(cost_rev_nl)) ) else - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, 0.0) + gen_cost[(n,i)] = JuMP.@expression(pm.model, 0.0) end end @@ -296,21 +297,21 @@ function _objective_min_fuel_and_flow_cost_polynomial_nl(pm::AbstractPowerModel; cost_rev = reverse(dcline["cost"]) if length(cost_rev) == 1 - dcline_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1]) + dcline_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1]) elseif length(cost_rev) == 2 - dcline_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1] + cost_rev[2]*p_dc) + dcline_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1] + cost_rev[2]*p_dc) elseif length(cost_rev) == 3 - dcline_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1] + cost_rev[2]*p_dc + cost_rev[3]*p_dc^2) + dcline_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1] + cost_rev[2]*p_dc + cost_rev[3]*p_dc^2) elseif length(cost_rev) >= 4 cost_rev_nl = cost_rev[4:end] - dcline_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1] + cost_rev[2]*p_dc + cost_rev[3]*p_dc^2 + sum( v*p_dc^(d+2) for (d,v) in enumerate(cost_rev_nl)) ) + dcline_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1] + cost_rev[2]*p_dc + cost_rev[3]*p_dc^2 + sum( v*p_dc^(d+2.0) for (d,v) in enumerate(cost_rev_nl)) ) else - dcline_cost[(n,i)] = JuMP.@NLexpression(pm.model, 0.0) + dcline_cost[(n,i)] = JuMP.@expression(pm.model, 0.0) end end end - return JuMP.@NLobjective(pm.model, Min, + return JuMP.@objective(pm.model, Min, sum( sum( gen_cost[(n,i)] for (i,gen) in nw_ref[:gen]) + sum( dcline_cost[(n,i)] for (i,dcline) in nw_ref[:dcline]) @@ -366,21 +367,22 @@ function _objective_min_fuel_cost_polynomial_nl(pm::AbstractPowerModel; report:: cost_rev = reverse(gen["cost"]) if length(cost_rev) == 1 - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1]) + gen_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1]) elseif length(cost_rev) == 2 - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1] + cost_rev[2]*pg) + gen_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1] + cost_rev[2]*pg) elseif length(cost_rev) == 3 - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1] + cost_rev[2]*pg + cost_rev[3]*pg^2) + gen_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1] + cost_rev[2]*pg + cost_rev[3]*pg^2) elseif length(cost_rev) >= 4 cost_rev_nl = cost_rev[4:end] - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, cost_rev[1] + cost_rev[2]*pg + cost_rev[3]*pg^2 + sum( v*pg^(d+3) for (d,v) in enumerate(cost_rev_nl)) ) + # TODO(odow): is this correct? Or should it be pg^(d+2.0)? + gen_cost[(n,i)] = JuMP.@expression(pm.model, cost_rev[1] + cost_rev[2]*pg + cost_rev[3]*pg^2 + sum( v*pg^(d+3.0) for (d,v) in enumerate(cost_rev_nl)) ) else - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, 0.0) + gen_cost[(n,i)] = JuMP.@expression(pm.model, 0.0) end end end - return JuMP.@NLobjective(pm.model, Min, + return JuMP.@objective(pm.model, Min, sum( sum( gen_cost[(n,i)] for (i,gen) in nw_ref[:gen] ) for (n, nw_ref) in nws(pm)) @@ -585,14 +587,14 @@ function objective_max_loadability(pm::AbstractPowerModel) time_elapsed = Dict(n => get(ref(pm, n), :time_elapsed, 1) for n in nws) load_weight = Dict(n => - Dict(i => get(load, "weight", 1.0) for (i,load) in ref(pm, n, :load)) + Dict(i => get(load, "weight", 1.0) for (i,load) in ref(pm, n, :load)) for n in nws) #println(load_weight) return JuMP.@objective(pm.model, Max, - sum( - ( + sum( + ( time_elapsed[n]*( sum(z_shunt[n][i] for (i,shunt) in ref(pm, n, :shunt)) + sum(load_weight[n][i]*abs(load["pd"])*z_demand[n][i] for (i,load) in ref(pm, n, :load)) diff --git a/src/form/acp.jl b/src/form/acp.jl index ddddf85a2..65dd99316 100644 --- a/src/form/acp.jl +++ b/src/form/acp.jl @@ -59,59 +59,27 @@ function constraint_power_balance(pm::AbstractACPModel, n::Int, i::Int, bus_arcs p_dc = get(var(pm, n), :p_dc, Dict()); _check_var_keys(p_dc, bus_arcs_dc, "active power", "dcline") q_dc = get(var(pm, n), :q_dc, Dict()); _check_var_keys(q_dc, bus_arcs_dc, "reactive power", "dcline") - # the check "typeof(p[arc]) <: JuMP.NonlinearExpression" is required for the - # case when p/q are nonlinear expressions instead of decision variables - # once NLExpressions are first order in JuMP it should be possible to - # remove this. - nl_form = length(bus_arcs) > 0 && (typeof(p[iterate(bus_arcs)[1]]) <: JuMP.NonlinearExpression) - - if !nl_form - cstr_p = JuMP.@constraint(pm.model, - sum(p[a] for a in bus_arcs) - + sum(p_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(psw[a_sw] for a_sw in bus_arcs_sw) - == - sum(pg[g] for g in bus_gens) - - sum(ps[s] for s in bus_storage) - - sum(pd for (i,pd) in bus_pd) - - sum(gs for (i,gs) in bus_gs)*vm^2 - ) - else - cstr_p = JuMP.@NLconstraint(pm.model, - sum(p[a] for a in bus_arcs) - + sum(p_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(psw[a_sw] for a_sw in bus_arcs_sw) - == - sum(pg[g] for g in bus_gens) - - sum(ps[s] for s in bus_storage) - - sum(pd for (i,pd) in bus_pd) - - sum(gs for (i,gs) in bus_gs)*vm^2 - ) - end + cstr_p = JuMP.@constraint(pm.model, + sum(p[a] for a in bus_arcs) + + sum(p_dc[a_dc] for a_dc in bus_arcs_dc) + + sum(psw[a_sw] for a_sw in bus_arcs_sw) + == + sum(pg[g] for g in bus_gens) + - sum(ps[s] for s in bus_storage) + - sum(pd for (i,pd) in bus_pd) + - sum(gs for (i,gs) in bus_gs)*vm^2 + ) - if !nl_form - cstr_q = JuMP.@constraint(pm.model, - sum(q[a] for a in bus_arcs) - + sum(q_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(qsw[a_sw] for a_sw in bus_arcs_sw) - == - sum(qg[g] for g in bus_gens) - - sum(qs[s] for s in bus_storage) - - sum(qd for (i,qd) in bus_qd) - + sum(bs for (i,bs) in bus_bs)*vm^2 - ) - else - cstr_q = JuMP.@NLconstraint(pm.model, - sum(q[a] for a in bus_arcs) - + sum(q_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(qsw[a_sw] for a_sw in bus_arcs_sw) - == - sum(qg[g] for g in bus_gens) - - sum(qs[s] for s in bus_storage) - - sum(qd for (i,qd) in bus_qd) - + sum(bs for (i,bs) in bus_bs)*vm^2 - ) - end + cstr_q = JuMP.@constraint(pm.model, + sum(q[a] for a in bus_arcs) + + sum(q_dc[a_dc] for a_dc in bus_arcs_dc) + + sum(qsw[a_sw] for a_sw in bus_arcs_sw) + == + sum(qg[g] for g in bus_gens) + - sum(qs[s] for s in bus_storage) + - sum(qd for (i,qd) in bus_qd) + + sum(bs for (i,bs) in bus_bs)*vm^2 + ) if _IM.report_duals(pm) sol(pm, n, :bus, i)[:lam_kcl_r] = cstr_p @@ -137,49 +105,26 @@ function constraint_power_balance_ls(pm::AbstractACPModel, n::Int, i::Int, bus_a z_shunt = get(var(pm, n), :z_shunt, Dict()); _check_var_keys(z_shunt, keys(bus_gs), "power factor", "shunt") # this is required for improved performance in NLP models - if length(z_shunt) <= 0 - cstr_p = JuMP.@constraint(pm.model, - sum(p[a] for a in bus_arcs) - + sum(p_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(psw[a_sw] for a_sw in bus_arcs_sw) - == - sum(pg[g] for g in bus_gens) - - sum(ps[s] for s in bus_storage) - - sum(pd*z_demand[i] for (i,pd) in bus_pd) - - sum(gs*z_shunt[i] for (i,gs) in bus_gs)*vm^2 - ) - cstr_q = JuMP.@constraint(pm.model, - sum(q[a] for a in bus_arcs) - + sum(q_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(qsw[a_sw] for a_sw in bus_arcs_sw) - == - sum(qg[g] for g in bus_gens) - - sum(qs[s] for s in bus_storage) - - sum(qd*z_demand[i] for (i,qd) in bus_qd) - + sum(bs*z_shunt[i] for (i,bs) in bus_bs)*vm^2 - ) - else - cstr_p = JuMP.@NLconstraint(pm.model, - sum(p[a] for a in bus_arcs) - + sum(p_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(psw[a_sw] for a_sw in bus_arcs_sw) - == - sum(pg[g] for g in bus_gens) - - sum(ps[s] for s in bus_storage) - - sum(pd*z_demand[i] for (i,pd) in bus_pd) - - sum(gs*z_shunt[i] for (i,gs) in bus_gs)*vm^2 - ) - cstr_q = JuMP.@NLconstraint(pm.model, - sum(q[a] for a in bus_arcs) - + sum(q_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(qsw[a_sw] for a_sw in bus_arcs_sw) - == - sum(qg[g] for g in bus_gens) - - sum(qs[s] for s in bus_storage) - - sum(qd*z_demand[i] for (i,qd) in bus_qd) - + sum(bs*z_shunt[i] for (i,bs) in bus_bs)*vm^2 - ) - end + cstr_p = JuMP.@constraint(pm.model, + sum(p[a] for a in bus_arcs) + + sum(p_dc[a_dc] for a_dc in bus_arcs_dc) + + sum(psw[a_sw] for a_sw in bus_arcs_sw) + == + sum(pg[g] for g in bus_gens) + - sum(ps[s] for s in bus_storage) + - sum(pd*z_demand[i] for (i,pd) in bus_pd) + - sum(gs*z_shunt[i] for (i,gs) in bus_gs)*vm^2 + ) + cstr_q = JuMP.@constraint(pm.model, + sum(q[a] for a in bus_arcs) + + sum(q_dc[a_dc] for a_dc in bus_arcs_dc) + + sum(qsw[a_sw] for a_sw in bus_arcs_sw) + == + sum(qg[g] for g in bus_gens) + - sum(qs[s] for s in bus_storage) + - sum(qd*z_demand[i] for (i,qd) in bus_qd) + + sum(bs*z_shunt[i] for (i,bs) in bus_bs)*vm^2 + ) if _IM.report_duals(pm) sol(pm, n, :bus, i)[:lam_kcl_r] = cstr_p @@ -237,8 +182,8 @@ function expression_branch_power_ohms_yt_from(pm::AbstractACPModel, n::Int, f_bu va_fr = var(pm, n, :va, f_bus) va_to = var(pm, n, :va, t_bus) - var(pm, n, :p)[f_idx] = JuMP.@NLexpression(pm.model, (g+g_fr)/tm^2*vm_fr^2 + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to)) ) - var(pm, n, :q)[f_idx] = JuMP.@NLexpression(pm.model, -(b+b_fr)/tm^2*vm_fr^2 - (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to)) ) + var(pm, n, :p)[f_idx] = JuMP.@expression(pm.model, (g+g_fr)/tm^2*vm_fr^2 + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to)) ) + var(pm, n, :q)[f_idx] = JuMP.@expression(pm.model, -(b+b_fr)/tm^2*vm_fr^2 - (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to)) ) end "" @@ -248,8 +193,8 @@ function expression_branch_power_ohms_yt_to(pm::AbstractACPModel, n::Int, f_bus, va_fr = var(pm, n, :va, f_bus) va_to = var(pm, n, :va, t_bus) - var(pm, n, :p)[t_idx] = JuMP.@NLexpression(pm.model, (g+g_to)*vm_to^2 + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr)) ) - var(pm, n, :q)[t_idx] = JuMP.@NLexpression(pm.model, -(b+b_to)*vm_to^2 - (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr)) ) + var(pm, n, :p)[t_idx] = JuMP.@expression(pm.model, (g+g_to)*vm_to^2 + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr)) ) + var(pm, n, :q)[t_idx] = JuMP.@expression(pm.model, -(b+b_to)*vm_to^2 - (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr)) ) end @@ -269,8 +214,8 @@ function constraint_ohms_yt_from(pm::AbstractACPModel, n::Int, f_bus, t_bus, f_i va_fr = var(pm, n, :va, f_bus) va_to = var(pm, n, :va, t_bus) - JuMP.@NLconstraint(pm.model, p_fr == (g+g_fr)/tm^2*vm_fr^2 + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to)) ) - JuMP.@NLconstraint(pm.model, q_fr == -(b+b_fr)/tm^2*vm_fr^2 - (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to)) ) + JuMP.@constraint(pm.model, p_fr == (g+g_fr)/tm^2*vm_fr^2 + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to)) ) + JuMP.@constraint(pm.model, q_fr == -(b+b_fr)/tm^2*vm_fr^2 - (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to)) ) end """ @@ -289,8 +234,8 @@ function constraint_ohms_yt_to(pm::AbstractACPModel, n::Int, f_bus, t_bus, f_idx va_fr = var(pm, n, :va, f_bus) va_to = var(pm, n, :va, t_bus) - JuMP.@NLconstraint(pm.model, p_to == (g+g_to)*vm_to^2 + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr)) ) - JuMP.@NLconstraint(pm.model, q_to == -(b+b_to)*vm_to^2 - (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr)) ) + JuMP.@constraint(pm.model, p_to == (g+g_to)*vm_to^2 + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr)) ) + JuMP.@constraint(pm.model, q_to == -(b+b_to)*vm_to^2 - (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr)) ) end """ @@ -309,8 +254,8 @@ function constraint_ohms_y_from(pm::AbstractACPModel, n::Int, f_bus, t_bus, f_id va_fr = var(pm, n, :va, f_bus) va_to = var(pm, n, :va, t_bus) - JuMP.@NLconstraint(pm.model, p_fr == (g+g_fr)*(vm_fr/tm)^2 - g*vm_fr/tm*vm_to*cos(va_fr-va_to-ta) + -b*vm_fr/tm*vm_to*sin(va_fr-va_to-ta) ) - JuMP.@NLconstraint(pm.model, q_fr == -(b+b_fr)*(vm_fr/tm)^2 + b*vm_fr/tm*vm_to*cos(va_fr-va_to-ta) + -g*vm_fr/tm*vm_to*sin(va_fr-va_to-ta) ) + JuMP.@constraint(pm.model, p_fr == (g+g_fr)*(vm_fr/tm)^2 - g*vm_fr/tm*vm_to*cos(va_fr-va_to-ta) + -b*vm_fr/tm*vm_to*sin(va_fr-va_to-ta) ) + JuMP.@constraint(pm.model, q_fr == -(b+b_fr)*(vm_fr/tm)^2 + b*vm_fr/tm*vm_to*cos(va_fr-va_to-ta) + -g*vm_fr/tm*vm_to*sin(va_fr-va_to-ta) ) end """ @@ -329,8 +274,8 @@ function constraint_ohms_y_to(pm::AbstractACPModel, n::Int, f_bus, t_bus, f_idx, va_fr = var(pm, n, :va, f_bus) va_to = var(pm, n, :va, t_bus) - JuMP.@NLconstraint(pm.model, p_to == (g+g_to)*vm_to^2 - g*vm_to*vm_fr/tm*cos(va_to-va_fr+ta) + -b*vm_to*vm_fr/tm*sin(va_to-va_fr+ta) ) - JuMP.@NLconstraint(pm.model, q_to == -(b+b_to)*vm_to^2 + b*vm_to*vm_fr/tm*cos(va_to-va_fr+ta) + -g*vm_to*vm_fr/tm*sin(va_to-va_fr+ta) ) + JuMP.@constraint(pm.model, p_to == (g+g_to)*vm_to^2 - g*vm_to*vm_fr/tm*cos(va_to-va_fr+ta) + -b*vm_to*vm_fr/tm*sin(va_to-va_fr+ta) ) + JuMP.@constraint(pm.model, q_to == -(b+b_to)*vm_to^2 + b*vm_to*vm_fr/tm*cos(va_to-va_fr+ta) + -g*vm_to*vm_fr/tm*sin(va_to-va_fr+ta) ) end @@ -382,8 +327,8 @@ function constraint_ohms_yt_from_on_off(pm::AbstractACPModel, n::Int, i, f_bus, va_to = var(pm, n, :va, t_bus) z = var(pm, n, :z_branch, i) - JuMP.@NLconstraint(pm.model, p_fr == z*( (g+g_fr)/tm^2*vm_fr^2 + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to))) ) - JuMP.@NLconstraint(pm.model, q_fr == z*(-(b+b_fr)/tm^2*vm_fr^2 - (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to))) ) + JuMP.@constraint(pm.model, p_fr == z*( (g+g_fr)/tm^2*vm_fr^2 + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to))) ) + JuMP.@constraint(pm.model, q_fr == z*(-(b+b_fr)/tm^2*vm_fr^2 - (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to))) ) end """ @@ -401,8 +346,8 @@ function constraint_ohms_yt_to_on_off(pm::AbstractACPModel, n::Int, i, f_bus, t_ va_to = var(pm, n, :va, t_bus) z = var(pm, n, :z_branch, i) - JuMP.@NLconstraint(pm.model, p_to == z*( (g+g_to)*vm_to^2 + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr))) ) - JuMP.@NLconstraint(pm.model, q_to == z*(-(b+b_to)*vm_to^2 - (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr))) ) + JuMP.@constraint(pm.model, p_to == z*( (g+g_to)*vm_to^2 + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr))) ) + JuMP.@constraint(pm.model, q_to == z*(-(b+b_to)*vm_to^2 - (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr))) ) end """ @@ -420,8 +365,8 @@ function constraint_ne_ohms_yt_from(pm::AbstractACPModel, n::Int, i, f_bus, t_bu va_to = var(pm, n, :va, t_bus) z = var(pm, n, :branch_ne, i) - JuMP.@NLconstraint(pm.model, p_fr == z*( (g+g_fr)/tm^2*vm_fr^2 + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to))) ) - JuMP.@NLconstraint(pm.model, q_fr == z*(-(b+b_fr)/tm^2*vm_fr^2 - (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to))) ) + JuMP.@constraint(pm.model, p_fr == z*( (g+g_fr)/tm^2*vm_fr^2 + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to))) ) + JuMP.@constraint(pm.model, q_fr == z*(-(b+b_fr)/tm^2*vm_fr^2 - (-b*tr-g*ti)/tm^2*(vm_fr*vm_to*cos(va_fr-va_to)) + (-g*tr+b*ti)/tm^2*(vm_fr*vm_to*sin(va_fr-va_to))) ) end """ @@ -439,8 +384,8 @@ function constraint_ne_ohms_yt_to(pm::AbstractACPModel, n::Int, i, f_bus, t_bus, va_to = var(pm, n, :va, t_bus) z = var(pm, n, :branch_ne, i) - JuMP.@NLconstraint(pm.model, p_to == z*( (g+g_to)*vm_to^2 + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr))) ) - JuMP.@NLconstraint(pm.model, q_to == z*(-(b+b_to)*vm_to^2 - (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr))) ) + JuMP.@constraint(pm.model, p_to == z*( (g+g_to)*vm_to^2 + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr))) ) + JuMP.@constraint(pm.model, q_to == z*(-(b+b_to)*vm_to^2 - (-b*tr+g*ti)/tm^2*(vm_to*vm_fr*cos(va_to-va_fr)) + (-g*tr-b*ti)/tm^2*(vm_to*vm_fr*sin(va_to-va_fr))) ) end """ @@ -460,8 +405,8 @@ function constraint_ohms_y_oltc_pst_from(pm::AbstractACPModel, n::Int, f_bus, t_ tm = var(pm, n, :tm, f_idx[1]) ta = var(pm, n, :ta, f_idx[1]) - JuMP.@NLconstraint(pm.model, p_fr == (g+g_fr)/tm^2*vm_fr^2 + (-g)/tm*(vm_fr*vm_to*cos(va_fr-va_to-ta)) + (-b)/tm*(vm_fr*vm_to*sin(va_fr-va_to-ta)) ) - JuMP.@NLconstraint(pm.model, q_fr == -(b+b_fr)/tm^2*vm_fr^2 - (-b)/tm*(vm_fr*vm_to*cos(va_fr-va_to-ta)) + (-g)/tm*(vm_fr*vm_to*sin(va_fr-va_to-ta)) ) + JuMP.@constraint(pm.model, p_fr == (g+g_fr)/tm^2*vm_fr^2 + (-g)/tm*(vm_fr*vm_to*cos(va_fr-va_to-ta)) + (-b)/tm*(vm_fr*vm_to*sin(va_fr-va_to-ta)) ) + JuMP.@constraint(pm.model, q_fr == -(b+b_fr)/tm^2*vm_fr^2 - (-b)/tm*(vm_fr*vm_to*cos(va_fr-va_to-ta)) + (-g)/tm*(vm_fr*vm_to*sin(va_fr-va_to-ta)) ) end """ @@ -481,8 +426,8 @@ function constraint_ohms_y_oltc_pst_to(pm::AbstractACPModel, n::Int, f_bus, t_bu tm = var(pm, n, :tm, f_idx[1]) ta = var(pm, n, :ta, f_idx[1]) - JuMP.@NLconstraint(pm.model, p_to == (g+g_to)*vm_to^2 + -g/tm*(vm_to*vm_fr*cos(va_to-va_fr+ta)) + -b/tm*(vm_to*vm_fr*sin(va_to-va_fr+ta)) ) - JuMP.@NLconstraint(pm.model, q_to == -(b+b_to)*vm_to^2 - -b/tm*(vm_to*vm_fr*cos(va_to-va_fr+ta)) + -g/tm*(vm_to*vm_fr*sin(va_to-va_fr+ta)) ) + JuMP.@constraint(pm.model, p_to == (g+g_to)*vm_to^2 + -g/tm*(vm_to*vm_fr*cos(va_to-va_fr+ta)) + -b/tm*(vm_to*vm_fr*sin(va_to-va_fr+ta)) ) + JuMP.@constraint(pm.model, q_to == -(b+b_to)*vm_to^2 - -b/tm*(vm_to*vm_fr*cos(va_to-va_fr+ta)) + -g/tm*(vm_to*vm_fr*sin(va_to-va_fr+ta)) ) end "`angmin <= z_branch[i]*(t[f_bus] - t[t_bus]) <= angmax`" @@ -549,13 +494,13 @@ function constraint_storage_losses(pm::AbstractACPModel, n::Int, i, bus, r, x, p sd = var(pm, n, :sd, i) qsc = var(pm, n, :qsc, i) - JuMP.@NLconstraint(pm.model, + JuMP.@constraint(pm.model, sum(ps[c] for c in conductors) + (sd - sc) == p_loss + sum(r[c]*(ps[c]^2 + qs[c]^2)/vm[c]^2 for c in conductors) ) - JuMP.@NLconstraint(pm.model, + JuMP.@constraint(pm.model, sum(qs[c] for c in conductors) == qsc + q_loss + sum(x[c]*(ps[c]^2 + qs[c]^2)/vm[c]^2 for c in conductors) diff --git a/src/form/acr.jl b/src/form/acr.jl index 2c6a93b10..61836dbcd 100644 --- a/src/form/acr.jl +++ b/src/form/acr.jl @@ -109,50 +109,26 @@ function constraint_power_balance_ls(pm::AbstractACRModel, n::Int, i::Int, bus_a z_demand = get(var(pm, n), :z_demand, Dict()); _check_var_keys(z_demand, keys(bus_pd), "power factor", "load") z_shunt = get(var(pm, n), :z_shunt, Dict()); _check_var_keys(z_shunt, keys(bus_gs), "power factor", "shunt") - # this is required for improved performance in NLP models - if length(z_shunt) <= 0 - cstr_p = JuMP.@constraint(pm.model, - sum(p[a] for a in bus_arcs) - + sum(p_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(psw[a_sw] for a_sw in bus_arcs_sw) - == - sum(pg[g] for g in bus_gens) - - sum(ps[s] for s in bus_storage) - - sum(pd*z_demand[i] for (i,pd) in bus_pd) - - sum(gs*z_shunt[i] for (i,gs) in bus_gs)*(vr^2 + vi^2) - ) - cstr_q = JuMP.@constraint(pm.model, - sum(q[a] for a in bus_arcs) - + sum(q_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(qsw[a_sw] for a_sw in bus_arcs_sw) - == - sum(qg[g] for g in bus_gens) - - sum(qs[s] for s in bus_storage) - - sum(qd*z_demand[i] for (i,qd) in bus_qd) - + sum(bs*z_shunt[i] for (i,bs) in bus_bs)*(vr^2 + vi^2) - ) - else - cstr_p = JuMP.@NLconstraint(pm.model, - sum(p[a] for a in bus_arcs) - + sum(p_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(psw[a_sw] for a_sw in bus_arcs_sw) - == - sum(pg[g] for g in bus_gens) - - sum(ps[s] for s in bus_storage) - - sum(pd*z_demand[i] for (i,pd) in bus_pd) - - sum(gs*z_shunt[i] for (i,gs) in bus_gs)*(vr^2 + vi^2) - ) - cstr_q = JuMP.@NLconstraint(pm.model, - sum(q[a] for a in bus_arcs) - + sum(q_dc[a_dc] for a_dc in bus_arcs_dc) - + sum(qsw[a_sw] for a_sw in bus_arcs_sw) - == - sum(qg[g] for g in bus_gens) - - sum(qs[s] for s in bus_storage) - - sum(qd*z_demand[i] for (i,qd) in bus_qd) - + sum(bs*z_shunt[i] for (i,bs) in bus_bs)*(vr^2 + vi^2) - ) - end + cstr_p = JuMP.@constraint(pm.model, + sum(p[a] for a in bus_arcs) + + sum(p_dc[a_dc] for a_dc in bus_arcs_dc) + + sum(psw[a_sw] for a_sw in bus_arcs_sw) + == + sum(pg[g] for g in bus_gens) + - sum(ps[s] for s in bus_storage) + - sum(pd*z_demand[i] for (i,pd) in bus_pd) + - sum(gs*z_shunt[i] for (i,gs) in bus_gs)*(vr^2 + vi^2) + ) + cstr_q = JuMP.@constraint(pm.model, + sum(q[a] for a in bus_arcs) + + sum(q_dc[a_dc] for a_dc in bus_arcs_dc) + + sum(qsw[a_sw] for a_sw in bus_arcs_sw) + == + sum(qg[g] for g in bus_gens) + - sum(qs[s] for s in bus_storage) + - sum(qd*z_demand[i] for (i,qd) in bus_qd) + + sum(bs*z_shunt[i] for (i,bs) in bus_bs)*(vr^2 + vi^2) + ) if _IM.report_duals(pm) sol(pm, n, :bus, i)[:lam_kcl_r] = cstr_p @@ -227,13 +203,13 @@ function constraint_storage_losses(pm::AbstractACRModel, n::Int, i, bus, r, x, p sd = var(pm, n, :sd, i) qsc = var(pm, n, :qsc, i) - JuMP.@NLconstraint(pm.model, + JuMP.@constraint(pm.model, sum(ps[c] for c in conductors) + (sd - sc) == p_loss + sum(r[c]*(ps[c]^2 + qs[c]^2)/(vr[c]^2 + vi[c]^2) for c in conductors) ) - JuMP.@NLconstraint(pm.model, + JuMP.@constraint(pm.model, sum(qs[c] for c in conductors) == qsc + q_loss + sum(x[c]*(ps[c]^2 + qs[c]^2)/(vr[c]^2 + vi[c]^2) for c in conductors) diff --git a/src/form/act.jl b/src/form/act.jl index 5124c2e92..7f5749e69 100644 --- a/src/form/act.jl +++ b/src/form/act.jl @@ -22,7 +22,7 @@ function constraint_model_voltage(pm::AbstractACTModel, n::Int) for (i,j) in ids(pm, n, :buspairs) JuMP.@constraint(pm.model, wr[(i,j)]^2 + wi[(i,j)]^2 == w[i]*w[j]) - JuMP.@NLconstraint(pm.model, wi[(i,j)] == tan(t[i] - t[j])*wr[(i,j)]) + JuMP.@constraint(pm.model, wi[(i,j)] == tan(t[i] - t[j])*wr[(i,j)]) end end diff --git a/src/form/iv.jl b/src/form/iv.jl index ca49d5bf2..4b820e7f9 100644 --- a/src/form/iv.jl +++ b/src/form/iv.jl @@ -50,8 +50,8 @@ function variable_gen_current(pm::AbstractIVRModel; nw::Int=nw_id_default, bound vi = var(pm, nw, :vi, busid) crg = var(pm, nw, :crg, i) cig = var(pm, nw, :cig, i) - pg[i] = JuMP.@NLexpression(pm.model, vr*crg + vi*cig) - qg[i] = JuMP.@NLexpression(pm.model, vi*crg - vr*cig) + pg[i] = JuMP.@expression(pm.model, vr*crg + vi*cig) + qg[i] = JuMP.@expression(pm.model, vi*crg - vr*cig) end var(pm, nw)[:pg] = pg var(pm, nw)[:qg] = qg @@ -85,10 +85,10 @@ function variable_dcline_current(pm::AbstractIVRModel; nw::Int=nw_id_default, bo cr_to = var(pm, nw, :crdc, (l,j,i)) ci_to = var(pm, nw, :cidc, (l,j,i)) - p[(l,i,j)] = JuMP.@NLexpression(pm.model, vr_fr*cr_fr + vi_fr*ci_fr) - q[(l,i,j)] = JuMP.@NLexpression(pm.model, vi_fr*cr_fr - vr_fr*ci_fr) - p[(l,j,i)] = JuMP.@NLexpression(pm.model, vr_to*cr_to + vi_to*ci_to) - q[(l,j,i)] = JuMP.@NLexpression(pm.model, vi_to*cr_to - vr_to*ci_to) + p[(l,i,j)] = JuMP.@expression(pm.model, vr_fr*cr_fr + vi_fr*ci_fr) + q[(l,i,j)] = JuMP.@expression(pm.model, vi_fr*cr_fr - vr_fr*ci_fr) + p[(l,j,i)] = JuMP.@expression(pm.model, vr_to*cr_to + vi_to*ci_to) + q[(l,j,i)] = JuMP.@expression(pm.model, vi_to*cr_to - vr_to*ci_to) end var(pm, nw)[:p_dc] = p @@ -189,14 +189,14 @@ function constraint_current_balance(pm::AbstractIVRModel, n::Int, i, bus_arcs, b crg = var(pm, n, :crg) cig = var(pm, n, :cig) - JuMP.@NLconstraint(pm.model, sum(cr[a] for a in bus_arcs) + JuMP.@constraint(pm.model, sum(cr[a] for a in bus_arcs) + sum(crdc[d] for d in bus_arcs_dc) == sum(crg[g] for g in bus_gens) - (sum(pd for pd in values(bus_pd))*vr + sum(qd for qd in values(bus_qd))*vi)/(vr^2 + vi^2) - sum(gs for gs in values(bus_gs))*vr + sum(bs for bs in values(bus_bs))*vi ) - JuMP.@NLconstraint(pm.model, sum(ci[a] for a in bus_arcs) + JuMP.@constraint(pm.model, sum(ci[a] for a in bus_arcs) + sum(cidc[d] for d in bus_arcs_dc) == sum(cig[g] for g in bus_gens) @@ -214,7 +214,7 @@ function constraint_thermal_limit_from(pm::AbstractIVRModel, n::Int, f_idx, rate crf = var(pm, n, :cr, f_idx) cif = var(pm, n, :ci, f_idx) - JuMP.@NLconstraint(pm.model, (vr^2 + vi^2)*(crf^2 + cif^2) <= rate_a^2) + JuMP.@constraint(pm.model, (vr^2 + vi^2)*(crf^2 + cif^2) <= rate_a^2) end "`p[t_idx]^2 + q[t_idx]^2 <= rate_a^2`" @@ -226,7 +226,7 @@ function constraint_thermal_limit_to(pm::AbstractIVRModel, n::Int, t_idx, rate_a crt = var(pm, n, :cr, t_idx) cit = var(pm, n, :ci, t_idx) - JuMP.@NLconstraint(pm.model, (vr^2 + vi^2)*(crt^2 + cit^2) <= rate_a^2) + JuMP.@constraint(pm.model, (vr^2 + vi^2)*(crt^2 + cit^2) <= rate_a^2) end """ @@ -386,15 +386,15 @@ function _objective_min_fuel_and_flow_cost_polynomial_linquad(pm::AbstractIVRMod for (i,gen) in nw_ref[:gen] bus = gen["gen_bus"] - #to avoid function calls inside of @NLconstraint: + #to avoid function calls inside of @constraint: pg = [var(pm, n, :pg, i) for c in conductor_ids(pm, n)] nc = length(conductor_ids(pm, n)) if length(gen["cost"]) == 1 gen_cost[(n,i)] = gen["cost"][1] elseif length(gen["cost"]) == 2 - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, gen["cost"][1]*sum(pg[c] for c in 1:nc) + gen["cost"][2]) + gen_cost[(n,i)] = JuMP.@expression(pm.model, gen["cost"][1]*sum(pg[c] for c in 1:nc) + gen["cost"][2]) elseif length(gen["cost"]) == 3 - gen_cost[(n,i)] = JuMP.@NLexpression(pm.model, gen["cost"][1]*sum(pg[c] for c in 1:nc)^2 + gen["cost"][2]*sum(pg[c] for c in 1:nc) + gen["cost"][3]) + gen_cost[(n,i)] = JuMP.@expression(pm.model, gen["cost"][1]*sum(pg[c] for c in 1:nc)^2 + gen["cost"][2]*sum(pg[c] for c in 1:nc) + gen["cost"][3]) else gen_cost[(n,i)] = 0.0 end @@ -403,23 +403,23 @@ function _objective_min_fuel_and_flow_cost_polynomial_linquad(pm::AbstractIVRMod from_idx = Dict(arc[1] => arc for arc in nw_ref[:arcs_from_dc]) for (i,dcline) in nw_ref[:dcline] bus = dcline["f_bus"] - #to avoid function calls inside of @NLconstraint: + #to avoid function calls inside of @constraint: p_dc = [var(pm, n, :p_dc, from_idx[i]) for c in conductor_ids(pm, n)] nc = length(conductor_ids(pm, n)) if length(dcline["cost"]) == 1 dcline_cost[(n,i)] = dcline["cost"][1] elseif length(dcline["cost"]) == 2 - dcline_cost[(n,i)] = JuMP.@NLexpression(pm.model, dcline["cost"][1]*sum(p_dc[c] for c in 1:nc) + dcline["cost"][2]) + dcline_cost[(n,i)] = JuMP.@expression(pm.model, dcline["cost"][1]*sum(p_dc[c] for c in 1:nc) + dcline["cost"][2]) elseif length(dcline["cost"]) == 3 - dcline_cost[(n,i)] = JuMP.@NLexpression(pm.model, dcline["cost"][1]*sum(p_dc[c] for c in 1:nc)^2 + dcline["cost"][2]*sum(p_dc[c] for c in 1:nc) + dcline["cost"][3]) + dcline_cost[(n,i)] = JuMP.@expression(pm.model, dcline["cost"][1]*sum(p_dc[c] for c in 1:nc)^2 + dcline["cost"][2]*sum(p_dc[c] for c in 1:nc) + dcline["cost"][3]) else dcline_cost[(n,i)] = 0.0 end end end - return JuMP.@NLobjective(pm.model, Min, + return JuMP.@objective(pm.model, Min, sum( sum( gen_cost[(n,i)] for (i,gen) in nw_ref[:gen] ) + sum( dcline_cost[(n,i)] for (i,dcline) in nw_ref[:dcline] ) @@ -433,7 +433,7 @@ function objective_variable_pg_cost(pm::AbstractIVRModel; report::Bool=true) for (n, nw_ref) in nws(pm) gen_lines = calc_cost_pwl_lines(nw_ref[:gen]) - #to avoid function calls inside of @NLconstraint + #to avoid function calls inside of @constraint pg_cost = var(pm, n)[:pg_cost] = JuMP.@variable(pm.model, [i in ids(pm, n, :gen)], base_name="$(n)_pg_cost", ) @@ -445,7 +445,7 @@ function objective_variable_pg_cost(pm::AbstractIVRModel; report::Bool=true) for (i, gen) in nw_ref[:gen] pg = [var(pm, n, :pg, i) for c in conductor_ids(pm, n)] for line in gen_lines[i] - JuMP.@NLconstraint(pm.model, pg_cost[i] >= line.slope*sum(pg[c] for c in 1:nc) + line.intercept) + JuMP.@constraint(pm.model, pg_cost[i] >= line.slope*sum(pg[c] for c in 1:nc) + line.intercept) end end end @@ -462,15 +462,15 @@ function objective_variable_dc_cost(pm::AbstractIVRModel; report::Bool=true) ) report && sol_component_value(pm, n, :dcline, :p_dc_cost, ids(pm, n, :dcline), dc_p_cost) - #to avoid function calls inside of @NLconstraint: + #to avoid function calls inside of @constraint: nc = length(conductor_ids(pm, n)) # dcline pwl cost for (i, dcline) in nw_ref[:dcline] arc = (i, dcline["f_bus"], dcline["t_bus"]) for line in dcline_lines[i] - #to avoid function calls inside of @NLconstraint: + #to avoid function calls inside of @constraint: p_dc = [var(pm, n, :p_dc, arc) for c in conductor_ids(pm, n)] - JuMP.@NLconstraint(pm.model, dc_p_cost[i] >= line.slope*sum(p_dc[c] for c in 1:nc) + line.intercept) + JuMP.@constraint(pm.model, dc_p_cost[i] >= line.slope*sum(p_dc[c] for c in 1:nc) + line.intercept) end end end diff --git a/test/docs.jl b/test/docs.jl index 9e4c1810b..e7b6eb891 100644 --- a/test/docs.jl +++ b/test/docs.jl @@ -30,7 +30,7 @@ #pretty print the model to the terminal #print(pm.model) - @test num_nonlinear_constraints(pm.model) == 12 + @test JuMP.num_nonlinear_constraints(pm.model) == 0 @test JuMP.num_variables(pm.model) == 28 result = optimize_model!(pm, optimizer=JuMP.optimizer_with_attributes(Ipopt.Optimizer, "print_level"=>0)) diff --git a/test/model.jl b/test/model.jl index 52b5eb330..c301e41c9 100644 --- a/test/model.jl +++ b/test/model.jl @@ -15,7 +15,7 @@ x = JuMP.@variable(m, my_var >= 0, start=0.0) pm = instantiate_model("../test/data/matpower/case5.m", ACPPowerModel, PowerModels.build_opf, jump_model=m) - @test num_nonlinear_constraints(pm.model) == 28 + @test JuMP.num_nonlinear_constraints(pm.model) == 0 @test JuMP.num_variables(pm.model) == 49 @test pm.model[:my_var] == x diff --git a/test/opf-var.jl b/test/opf-var.jl index e6ee6b862..4908933c5 100644 --- a/test/opf-var.jl +++ b/test/opf-var.jl @@ -127,8 +127,8 @@ end result = PowerModels._solve_opf_cl(data, SDPWRMPowerModel, sdp_solver) @test result["termination_status"] == OPTIMAL - @test isapprox(result["objective"], 5728.62; atol = 1e0) - #@test isapprox(result["objective"], 5747.63; atol = 1e0) + # @test isapprox(result["objective"], 5728.62; atol = 1e0) + @test isapprox(result["objective"], 5747.63; atol = 1e0) end @testset "5-bus case" begin data = build_current_data("../test/data/matpower/case5.m") @@ -144,8 +144,8 @@ end result = PowerModels._solve_opf_cl(data, SDPWRMPowerModel, sdp_solver) @test result["termination_status"] == OPTIMAL || result["termination_status"] == ALMOST_OPTIMAL - @test isapprox(result["objective"], 7505.33; atol = 1e0) - #@test isapprox(result["objective"], 7637.95; atol = 1e0) + # @test isapprox(result["objective"], 7505.33; atol = 1e0) + @test isapprox(result["objective"], 7637.95; atol = 1e0) end end diff --git a/test/opf.jl b/test/opf.jl index 5f8419ba5..36a860fa4 100644 --- a/test/opf.jl +++ b/test/opf.jl @@ -590,8 +590,8 @@ end result = run_opf("../test/data/matpower/case3.m", SOCWRConicPowerModel, sdp_solver) @test result["termination_status"] == OPTIMAL - @test isapprox(result["objective"], 5736.94; atol = 2e0) - #@test isapprox(result["objective"], 5747.37; atol = 2e0) + # @test isapprox(result["objective"], 5736.94; atol = 2e0) + @test isapprox(result["objective"], 5747.37; atol = 2e0) end @testset "5-bus transformer swap case" begin result = run_opf("../test/data/matpower/case5.m", SOCWRConicPowerModel, sdp_solver) @@ -621,31 +621,31 @@ end result = run_opf("../test/data/matpower/case5_npg.m", SOCWRConicPowerModel, sdp_solver) @test result["termination_status"] == OPTIMAL - @test isapprox(result["objective"], 3551.71; atol = 40) - #@test isapprox(result["objective"], 3602.11; atol = 40) + # @test isapprox(result["objective"], 3551.71; atol = 40) + @test isapprox(result["objective"], 3602.11; atol = 40) end @testset "5-bus with pwl costs" begin result = run_opf("../test/data/matpower/case5_pwlc.m", SOCWRConicPowerModel, sdp_solver) @test result["termination_status"] == OPTIMAL - @test isapprox(result["objective"], 42889; atol = 1e0) - #@test isapprox(result["objective"], 42906; atol = 1e0) + # @test isapprox(result["objective"], 42889; atol = 1e0) + @test isapprox(result["objective"], 42906; atol = 1e0) end @testset "6-bus case" begin result = run_opf("../test/data/matpower/case6.m", SOCWRConicPowerModel, sdp_solver) @test result["termination_status"] == OPTIMAL #@test isapprox(result["objective"], 11472.2; atol = 3e0) - @test isapprox(result["objective"], 11451.5; atol = 3e0) - #@test isapprox(result["objective"], 11473.4; atol = 3e0) + # @test isapprox(result["objective"], 11451.5; atol = 3e0) + @test isapprox(result["objective"], 11473.4; atol = 3e0) end @testset "24-bus rts case" begin result = run_opf("../test/data/matpower/case24.m", SOCWRConicPowerModel, sdp_solver) @test result["termination_status"] == OPTIMAL #@test isapprox(result["objective"], 70693.9; atol = 1e0) - @test isapprox(result["objective"], 70670.0; atol = 1e0) - #@test isapprox(result["objective"], 70683.5; atol = 1e0) + # @test isapprox(result["objective"], 70670.0; atol = 1e0) + @test isapprox(result["objective"], 70683.5; atol = 1e0) end @testset "14-bus variable bounds" begin pm = instantiate_model("../test/data/matpower/case14.m", SOCWRConicPowerModel, PowerModels.build_opf) @@ -906,14 +906,16 @@ end result = run_opf("../test/data/matpower/case3.m", SDPWRMPowerModel, sdp_solver) @test result["termination_status"] == OPTIMAL - @test isapprox(result["objective"], 5818.00; atol = 1e1) - #@test isapprox(result["objective"], 5852.51; atol = 1e1) + # @test isapprox(result["objective"], 5818.00; atol = 1e1) + @test isapprox(result["objective"], 5852.51; atol = 1e1) @test haskey(result["solution"],"WR") @test haskey(result["solution"],"WI") - @test isapprox(result["solution"]["bus"]["1"]["w"], 1.179, atol = 1e-2) + # @test isapprox(result["solution"]["bus"]["1"]["w"], 1.179, atol = 1e-2) + @test isapprox(result["solution"]["bus"]["1"]["w"], 1.210, atol = 1e-2) @test isapprox(result["solution"]["branch"]["1"]["wr"], 0.941, atol = 1e-2) - @test isapprox(result["solution"]["branch"]["1"]["wi"], 0.269, atol = 1e-2) + # @test isapprox(result["solution"]["branch"]["1"]["wi"], 0.269, atol = 1e-2) + @test isapprox(result["solution"]["branch"]["1"]["wi"], 0.285, atol = 1e-2) end @testset "5-bus asymmetric case" begin result = run_opf("../test/data/matpower/case5_asym.m", SDPWRMPowerModel, sdp_solver) @@ -938,8 +940,8 @@ end @test result["termination_status"] == OPTIMAL #@test isapprox(result["objective"], 6827.34; atol = 1e0) - @test isapprox(result["objective"], 6735.17; atol = 1e0) - #@test isapprox(result["objective"], 6827.71; atol = 1e0) + # @test isapprox(result["objective"], 6735.17; atol = 1e0) + @test isapprox(result["objective"], 6827.71; atol = 1e0) end # too slow for unit tests # @testset "14-bus case" begin @@ -953,8 +955,8 @@ end @test result["termination_status"] == OPTIMAL #@test isapprox(result["objective"], 11580.8; atol = 1e1) - @test isapprox(result["objective"], 11507.7; atol = 1e1) - #@test isapprox(result["objective"], 11580.5; atol = 1e1) + # @test isapprox(result["objective"], 11507.7; atol = 1e1) + @test isapprox(result["objective"], 11580.5; atol = 1e1) end @testset "14-bus variable bounds" begin pm = instantiate_model("../test/data/matpower/case14.m", SDPWRMPowerModel, PowerModels.build_opf) @@ -969,14 +971,16 @@ end @test result["termination_status"] == OPTIMAL #@test isapprox(result["objective"], 5851.23; atol = 1e1) - @test isapprox(result["objective"], 5818.00; atol = 1e1) - #@test isapprox(result["objective"], 5852.51; atol = 1e1) + # @test isapprox(result["objective"], 5818.00; atol = 1e1) + @test isapprox(result["objective"], 5852.51; atol = 1e1) @test haskey(result["solution"]["w_group"]["1"],"WR") @test haskey(result["solution"]["w_group"]["1"],"WI") - @test isapprox(result["solution"]["bus"]["1"]["w"], 1.179, atol = 1e-2) + # @test isapprox(result["solution"]["bus"]["1"]["w"], 1.179, atol = 1e-2) + @test isapprox(result["solution"]["bus"]["1"]["w"], 1.210, atol = 1e-2) @test isapprox(result["solution"]["branch"]["1"]["wr"], 0.941, atol = 1e-2) - @test isapprox(result["solution"]["branch"]["1"]["wi"], 0.269, atol = 1e-2) + # @test isapprox(result["solution"]["branch"]["1"]["wi"], 0.269, atol = 1e-2) + @test isapprox(result["solution"]["branch"]["1"]["wi"], 0.285, atol = 1e-2) end @testset "5-bus with asymmetric line charge" begin result = run_opf("../test/data/pti/case5_alc.raw", SparseSDPWRMPowerModel, sdp_solver) diff --git a/test/runtests.jl b/test/runtests.jl index 49c7902f0..f1f9f59f6 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,3 +1,6 @@ +import Pkg +Pkg.pkg"add JuMP#od/nlp-expr MathOptInterface#od/nlp-expr Ipopt#od/nlp-expr" + using PowerModels import InfrastructureModels import Memento @@ -18,16 +21,6 @@ import LinearAlgebra import SparseArrays using Test - -# compat for JuMP v0.22/v0.23 transition -# can be removed after dropping support for v0.22 -if !isdefined(JuMP, :num_nonlinear_constraints) - num_nonlinear_constraints = JuMP.num_nl_constraints -else - num_nonlinear_constraints = JuMP.num_nonlinear_constraints -end - - # default setup for solvers nlp_solver = JuMP.optimizer_with_attributes(Ipopt.Optimizer, "tol"=>1e-6, "print_level"=>0) nlp_ws_solver = JuMP.optimizer_with_attributes(Ipopt.Optimizer, "tol"=>1e-6, "mu_init"=>1e-4, "print_level"=>0)