From 27f2f1c1d4c873ce3bd78c7307bd2c9b2c482a40 Mon Sep 17 00:00:00 2001 From: Joaquim Dias Garcia Date: Mon, 3 Oct 2022 19:02:40 -0300 Subject: [PATCH 1/8] Fix sign in lower level maximization QP Close #182 Needs tests --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 456d9b74..6e3f60e4 100644 --- a/Project.toml +++ b/Project.toml @@ -10,7 +10,7 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" [compat] -Dualization = "0.5.4" +Dualization = "0.5.5" JuMP = "1.0" MathOptInterface = "1.2" julia = "1.6" From d3cb749a8581ee0920d881b8dead1e15326a182a Mon Sep 17 00:00:00 2001 From: Joaquim Dias Garcia Date: Mon, 3 Oct 2022 19:31:13 -0300 Subject: [PATCH 2/8] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 6e3f60e4..ebc05638 100644 --- a/Project.toml +++ b/Project.toml @@ -11,6 +11,6 @@ MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" [compat] Dualization = "0.5.5" -JuMP = "1.0" +JuMP = "1.0.0 - 1.1.1" MathOptInterface = "1.2" julia = "1.6" From 1b2a8f420ea8c62bafc411c2ae5f9d48419449e7 Mon Sep 17 00:00:00 2001 From: Joaquim Dias Garcia Date: Mon, 3 Oct 2022 20:08:42 -0300 Subject: [PATCH 3/8] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index ebc05638..1ba4e754 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "BilevelJuMP" uuid = "485130c0-026e-11ea-0f1a-6992cd14145c" authors = ["Joaquim Garcia , guilhermebodin "] -version = "0.5.1" +version = "0.5.2" [deps] Dualization = "191a621a-6537-11e9-281d-650236a99e60" From 1b48de7c23cb9cffcc3f41c487785d661e402f56 Mon Sep 17 00:00:00 2001 From: joaquim Date: Mon, 3 Oct 2022 20:22:56 -0300 Subject: [PATCH 4/8] clean deps --- docs/Project.toml | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/docs/Project.toml b/docs/Project.toml index ee60097f..5664137e 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -1,14 +1,12 @@ [deps] -BilevelJuMP = "485130c0-026e-11ea-0f1a-6992cd14145c" Cbc = "9961bab8-2fa3-5c5a-9d89-47fab24efd76" Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" Ipopt = "b6b21f68-93f8-5de0-b562-5493be1d77c9" -JuMP = "4076af6c-e467-56ae-b986-b466b2749572" Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306" MibS_jll = "118347d2-e127-56b9-9933-6abf9cc1adc1" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [compat] Documenter = "~0.26" -MibS_jll = "=1.1.3" Ipopt = "=1.0.2" +MibS_jll = "=1.1.3" From 13d96ef4cfb6c3d9cdb7229186394334c04b9e62 Mon Sep 17 00:00:00 2001 From: joaquim Date: Tue, 4 Oct 2022 02:13:28 -0300 Subject: [PATCH 5/8] fix doc project --- docs/Project.toml | 2 ++ 1 file changed, 2 insertions(+) diff --git a/docs/Project.toml b/docs/Project.toml index 5664137e..9757f6bf 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -1,7 +1,9 @@ [deps] +BilevelJuMP = "485130c0-026e-11ea-0f1a-6992cd14145c" Cbc = "9961bab8-2fa3-5c5a-9d89-47fab24efd76" Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" Ipopt = "b6b21f68-93f8-5de0-b562-5493be1d77c9" +JuMP = "4076af6c-e467-56ae-b986-b466b2749572" Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306" MibS_jll = "118347d2-e127-56b9-9933-6abf9cc1adc1" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" From 9344cf2121c80ff8b5f24667af75700bb8d9227c Mon Sep 17 00:00:00 2001 From: joaquim Date: Tue, 4 Oct 2022 02:37:34 -0300 Subject: [PATCH 6/8] add max ll maxqp tests --- test/jump.jl | 168 +++++++++++++++++++++++++++++++++++++---------- test/runtests.jl | 61 +++++++++-------- 2 files changed, 168 insertions(+), 61 deletions(-) diff --git a/test/jump.jl b/test/jump.jl index 77632220..9d7b0346 100644 --- a/test/jump.jl +++ b/test/jump.jl @@ -1488,7 +1488,7 @@ end # TODO - add quadratic problems # 9.3.2 - parg 221 -function jump_HTP_quad01(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Config()) +function jump_HTP_quad01(optimizer, is_min, mode = BilevelJuMP.SOS1Mode(), config = Config()) atol = config.atol @@ -1503,8 +1503,13 @@ function jump_HTP_quad01(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf @constraint(Upper(model), x >= 0) @constraint(Upper(model), y >= 0) # only in lowrrin GAMS - @objective(Lower(model), Min, - (y-1)^2 -1.5*x*y) + if is_min + @objective(Lower(model), Min, + (y-1)^2 -1.5*x*y) + else + @objective(Lower(model), Max, + -((y-1)^2 -1.5*x*y)) + end @constraint(Lower(model), -3x + y <= -3) @constraint(Lower(model), x - 0.5y <= 4) @@ -1523,7 +1528,7 @@ function jump_HTP_quad01(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf end # 9.3.3- parg 222 -function jump_HTP_quad02(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Config()) +function jump_HTP_quad02(optimizer, is_min, mode = BilevelJuMP.SOS1Mode(), config = Config()) atol = config.atol @@ -1541,8 +1546,13 @@ function jump_HTP_quad02(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf @constraint(Upper(model), y >= 0) @constraint(Upper(model), y <= 20) - @objective(Lower(model), Min, - (x + 2y - 30)^2) + if is_min + @objective(Lower(model), Min, + (x + 2y - 30)^2) + else + @objective(Lower(model), Max, + -((x + 2y - 30)^2)) + end @constraint(Lower(model), x + y <= 20) @constraint(Lower(model), - y <= 0) @@ -1560,7 +1570,7 @@ function jump_HTP_quad02(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf end # 9.3.4- parg 223 -function jump_HTP_quad03(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Config()) +function jump_HTP_quad03(optimizer, is_min, mode = BilevelJuMP.SOS1Mode(), config = Config()) atol = config.atol start = config.start_value @@ -1579,12 +1589,17 @@ function jump_HTP_quad03(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf @constraint(Upper(model), [i=1:2], y[i] >= -10) @constraint(Upper(model), [i=1:2], y[i] <= 20) - @objective(Lower(model), Min, - (-x[1] + y[1] + 40)^2 + (-x[2] + y[2] + 20)^2) - # the boo does not contain the 40, it is a 20 there - # however the solution does not match - # the file https://www.gams.com/latest/emplib_ml/libhtml/emplib_flds923.html - # has a 40 + if is_min + @objective(Lower(model), Min, + (-x[1] + y[1] + 40)^2 + (-x[2] + y[2] + 20)^2) + else + @objective(Lower(model), Max, + -((-x[1] + y[1] + 40)^2 + (-x[2] + y[2] + 20)^2)) + end + # the boo does not contain the 40, it is a 20 there + # however the solution does not match + # the file https://www.gams.com/latest/emplib_ml/libhtml/emplib_flds923.html + # has a 40 @constraint(Lower(model), [i=1:2],- x[i] + 2y[i] <= -10) @constraint(Lower(model), [i=1:2], y[i] >= -10) @@ -1604,7 +1619,7 @@ function jump_HTP_quad03(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf end # 9.3.5- parg 225 -function jump_HTP_quad04(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Config()) +function jump_HTP_quad04(optimizer, is_min, mode = BilevelJuMP.SOS1Mode(), config = Config()) atol = config.atol start = config.start_value @@ -1619,8 +1634,13 @@ function jump_HTP_quad04(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf 0.5*((y[1]-2)^2+(y[2]-2)^2) ) @constraint(Upper(model), [i=1:2], y[i] >= 0) - @objective(Lower(model), Min, - 0.5*(y[1]^2) + y[2]) + if is_min + @objective(Lower(model), Min, + 0.5*(y[1]^2) + y[2]) + else + @objective(Lower(model), Max, + -(0.5*(y[1]^2) + y[2])) + end @constraint(Lower(model), y[1] + y[2] == x) @constraint(Lower(model), [i=1:2], y[i] >= 0) @@ -1636,7 +1656,7 @@ function jump_HTP_quad04(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf end # 9.3.6- parg 226 -function jump_HTP_quad05(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Config()) +function jump_HTP_quad05(optimizer, is_min, mode = BilevelJuMP.SOS1Mode(), config = Config()) atol = config.atol start = config.start_value @@ -1652,8 +1672,13 @@ function jump_HTP_quad05(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf @constraint(Upper(model), x >= 0) @constraint(Upper(model), x <= 8) - @objective(Lower(model), Min, - (y-5)^2 ) + if is_min + @objective(Lower(model), Min, + (y-5)^2 ) + else + @objective(Lower(model), Max, + -((y-5)^2) ) + end @constraint(Lower(model), -2x+y <= 1) @constraint(Lower(model), x-2y <= 2) @@ -1672,7 +1697,7 @@ function jump_HTP_quad05(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf end # 9.3.7- parg 227 -function jump_HTP_quad06(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Config()) +function jump_HTP_quad06(optimizer, is_min, mode = BilevelJuMP.SOS1Mode(), config = Config()) atol = config.atol start = config.start_value @@ -1687,8 +1712,13 @@ function jump_HTP_quad06(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf x[1]^2 -2x[1] +x[2]^2 -2x[2] +y[1]^2 +y[2]^2) @constraint(Upper(model), [i=1:2], y[i] >= 0) - @objective(Lower(model), Min, - (-x[1] + y[1])^2 + (-x[2] + y[2])^2) + if is_min + @objective(Lower(model), Min, + (-x[1] + y[1])^2 + (-x[2] + y[2])^2) + else + @objective(Lower(model), Max, + -((-x[1] + y[1])^2 + (-x[2] + y[2])^2)) + end @constraint(Lower(model), [i=1:2], y[i] >= 0.5) @constraint(Lower(model), [i=1:2], y[i] <= 1.5) @@ -1704,7 +1734,7 @@ function jump_HTP_quad06(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf @test sol ≈ [0.5 ; 0.5; 0.5; 0.5] atol=atol end -function jump_HTP_quad06b(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Config()) +function jump_HTP_quad06b(optimizer, is_min, mode = BilevelJuMP.SOS1Mode(), config = Config()) # TODO reviews the behaviour comment atol = config.atol start = config.start_value @@ -1719,8 +1749,13 @@ function jump_HTP_quad06b(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Con x[1]^2 +x[2]^2 +y[1]^2 - 3y[1] +y[2]^2 - 3y[2]) @constraint(Upper(model), [i=1:2], y[i] >= 0) - @objective(Lower(model), Min, - (-x[1] + y[1])^2 + (-x[2] + y[2])^2) + if is_min + @objective(Lower(model), Min, + (-x[1] + y[1])^2 + (-x[2] + y[2])^2) + else + @objective(Lower(model), Max, + -((-x[1] + y[1])^2 + (-x[2] + y[2])^2)) + end @constraint(Lower(model), [i=1:2], y[i] >= 0.5) @constraint(Lower(model), [i=1:2], y[i] <= 1.5) @@ -1755,7 +1790,7 @@ end # 9.3.8- parg 228 -function jump_HTP_quad07(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Config()) +function jump_HTP_quad07(optimizer, is_min, mode = BilevelJuMP.SOS1Mode(), config = Config()) atol = config.atol @@ -1770,8 +1805,13 @@ function jump_HTP_quad07(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf @constraint(Upper(model), x >= 0) @constraint(Upper(model), y >= 0) # only in lowrrin GAMS - @objective(Lower(model), Min, - (y-1)^2 -1.5*x*y) + if is_min + @objective(Lower(model), Min, + (y-1)^2 -1.5*x*y) + else + @objective(Lower(model), Max, + -((y-1)^2 -1.5*x*y)) + end @constraint(Lower(model), -3x + y <= -3) @constraint(Lower(model), x - 0.5y <= 4) @@ -1790,7 +1830,7 @@ function jump_HTP_quad07(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf end # 9.3.9 - parg 229 -function jump_HTP_quad08(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Config()) +function jump_HTP_quad08(optimizer, is_min, mode = BilevelJuMP.SOS1Mode(), config = Config()) # Q objective is not PSD atol = config.atol @@ -1809,8 +1849,13 @@ function jump_HTP_quad08(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf @constraint(Upper(model), y >= 0) @constraint(Upper(model), y <= 1) - @objective(Lower(model), Min, - -(1-4x)*y -(2x+2) ) + if is_min + @objective(Lower(model), Min, + -(1-4x)*y -(2x+2) ) + else + @objective(Lower(model), Max, + -(-(1-4x)*y -(2x+2) )) + end @constraint(Lower(model), y >= 0) @constraint(Lower(model), y <= 1) @@ -1827,7 +1872,7 @@ function jump_HTP_quad08(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf end # 9.3.10- parg 230 -function jump_HTP_quad09(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Config()) +function jump_HTP_quad09(optimizer, is_min, mode = BilevelJuMP.SOS1Mode(), config = Config()) atol = config.atol start = config.start_value @@ -1843,8 +1888,13 @@ function jump_HTP_quad09(optimizer, mode = BilevelJuMP.SOS1Mode(), config = Conf @constraint(Upper(model), [i=1:2], y[i] >= 0) @constraint(Upper(model), x >= 0) - @objective(Lower(model), Min, - 2y[1]+x*y[2]) + if is_min + @objective(Lower(model), Min, + 2y[1]+x*y[2]) + else + @objective(Lower(model), Max, + -(2y[1]+x*y[2])) + end @constraint(Lower(model), x + 4 <= y[1] + y[2]) # this 4 is missing from the book @@ -2826,4 +2876,54 @@ function jump_01_sum_agg(optimizer, config = Config()) # @test dual(c3) ≈ [0] atol=atol # @test dual(c4) ≈ [1] atol=atol #NLP fail +end + +""" +From: https://github.com/joaquimg/BilevelJuMP.jl/issues/182 +By: LukasBarner + +Siddiqui S, Gabriel SA (2013). An sos1-based approach for solving mpecs with a natural gas market applica- +tion. Networks and Spatial Economics 13(2):205–227. +""" +function jump_qp_lower_min(config = Config()) + F = [1,2] + c = Dict(1=>1, 2=>1) + C = 1 + a = 13 + b = 1 + + model = BilevelModel(Ipopt.Optimizer, mode = BilevelJuMP.ProductMode()) + + @variable(Lower(model), q[F] >= 0) + @variable(Upper(model), Q >= 0) + + @objective(Upper(model), Max, ((a-b * (q[1] + q[2] + Q)) * Q - C*Q) ) + + @objective(Lower(model), Min, -((a-b * (q[1] + q[2] + Q)) * q[1] - C*q[1] + (a-b * (q[1] + q[2] + Q)) * q[2] - C*q[2] + b*q[1]*q[2]) ) + + optimize!(model) + + @test isapprox(value.(Q), 6; atol=1e-5) + @test isapprox(value.(q).data, [2,2]; atol=1e-5) +end +function jump_qp_lower_min(config = Config()) + F = [1,2] + c = Dict(1=>1, 2=>1) + C = 1 + a = 13 + b = 1 + + model = BilevelModel(Ipopt.Optimizer, mode = BilevelJuMP.ProductMode()) + + @variable(Lower(model), q[F] >= 0) + @variable(Upper(model), Q >= 0) + + @objective(Upper(model), Max, ((a-b * (q[1] + q[2] + Q)) * Q - C*Q) ) + + @objective(Lower(model), Max, +((a-b * (q[1] + q[2] + Q)) * q[1] - C*q[1] + (a-b * (q[1] + q[2] + Q)) * q[2] - C*q[2] + b*q[1]*q[2]) ) + + optimize!(model) + + @test isapprox(value.(Q), 6; atol=1e-5) + @test isapprox(value.(q).data, [2,2]; atol=1e-5) end \ No newline at end of file diff --git a/test/runtests.jl b/test/runtests.jl index fba42199..ca1ae9cb 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -280,34 +280,36 @@ end end @testset "Princeton Handbook Quadratic" begin - for solver in vcat(solvers_nlp, solvers_nlp_sum) - jump_HTP_quad01(solver.opt, solver.mode) - jump_HTP_quad02(solver.opt, solver.mode) - jump_HTP_quad04(solver.opt, solver.mode, CONFIG_3) - jump_HTP_quad05(solver.opt, solver.mode) - jump_HTP_quad06(solver.opt, solver.mode, CONFIG_3) - # jump_HTP_quad06b(solver.opt, solver.mode) - jump_HTP_quad07(solver.opt, solver.mode) - jump_HTP_quad08(solver.opt, solver.mode) # not PSD - jump_HTP_quad09(solver.opt, solver.mode) - end - for solver in solvers_nlp - jump_HTP_quad03(solver.opt, solver.mode) - end + for is_min in [true, false] + for solver in vcat(solvers_nlp, solvers_nlp_sum) + jump_HTP_quad01(solver.opt, is_min, solver.mode) + jump_HTP_quad02(solver.opt, is_min, solver.mode) + jump_HTP_quad04(solver.opt, is_min, solver.mode, CONFIG_3) + jump_HTP_quad05(solver.opt, is_min, solver.mode) + jump_HTP_quad06(solver.opt, is_min, solver.mode, CONFIG_3) + # jump_HTP_quad06b(solver.opt, is_min, solver.mode) + jump_HTP_quad07(solver.opt, is_min, solver.mode) + jump_HTP_quad08(solver.opt, is_min, solver.mode) # not PSD + jump_HTP_quad09(solver.opt, is_min, solver.mode) + end + for solver in solvers_nlp + jump_HTP_quad03(solver.opt, is_min, solver.mode) + end - for solver in solvers_sos_quad - jump_HTP_quad01(solver.opt, solver.mode, CONFIG_4) - jump_HTP_quad02(solver.opt, solver.mode) - jump_HTP_quad04(solver.opt, solver.mode) - jump_HTP_quad05(solver.opt, solver.mode) - jump_HTP_quad06(solver.opt, solver.mode) - jump_HTP_quad06b(solver.opt, solver.mode, CONFIG_4) - jump_HTP_quad07(solver.opt, solver.mode, CONFIG_4) - # jump_HTP_quad08(solver.opt, solver.mode) # not PSD - end - for solver in solvers_sos - jump_HTP_quad03(solver.opt, solver.mode) - jump_HTP_quad09(solver.opt, solver.mode) + for solver in solvers_sos_quad + jump_HTP_quad01(solver.opt, is_min, solver.mode, CONFIG_4) + jump_HTP_quad02(solver.opt, is_min, solver.mode) + jump_HTP_quad04(solver.opt, is_min, solver.mode) + jump_HTP_quad05(solver.opt, is_min, solver.mode) + jump_HTP_quad06(solver.opt, is_min, solver.mode) + jump_HTP_quad06b(solver.opt, is_min, solver.mode, CONFIG_4) + jump_HTP_quad07(solver.opt, is_min, solver.mode, CONFIG_4) + # jump_HTP_quad08(solver.opt, is_min, solver.mode) # not PSD + end + for solver in solvers_sos + jump_HTP_quad03(solver.opt, is_min, solver.mode) + jump_HTP_quad09(solver.opt, is_min, solver.mode) + end end end @@ -368,6 +370,11 @@ end end end +@testset "Lower QP" begin + jump_qp_lower_min() + jump_qp_lower_max() +end + @testset "Fruits" begin for solver in solvers_fa2 jump_fruits(solver.opt, solver.mode, CONFIG_4, 0.05) From a2a08b798a80e91cf912f55f11d47e4a8d0424da Mon Sep 17 00:00:00 2001 From: joaquim Date: Tue, 4 Oct 2022 08:33:59 -0300 Subject: [PATCH 7/8] fix failing test --- test/jump.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/jump.jl b/test/jump.jl index 9d7b0346..0e6afddd 100644 --- a/test/jump.jl +++ b/test/jump.jl @@ -2906,7 +2906,7 @@ function jump_qp_lower_min(config = Config()) @test isapprox(value.(Q), 6; atol=1e-5) @test isapprox(value.(q).data, [2,2]; atol=1e-5) end -function jump_qp_lower_min(config = Config()) +function jump_qp_lower_max(config = Config()) F = [1,2] c = Dict(1=>1, 2=>1) C = 1 From 03e8c0bcefa71d3870eb10bfdae63512c86af739 Mon Sep 17 00:00:00 2001 From: joaquim Date: Tue, 4 Oct 2022 10:59:02 -0300 Subject: [PATCH 8/8] bump dualization to fix parameters --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 1ba4e754..1d0bc014 100644 --- a/Project.toml +++ b/Project.toml @@ -10,7 +10,7 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" [compat] -Dualization = "0.5.5" +Dualization = "0.5.6" JuMP = "1.0.0 - 1.1.1" MathOptInterface = "1.2" julia = "1.6"