Skip to content

Commit dc17a6e

Browse files
dnguyen227dnguyen227
authored
Support cutting planes technique. (#125)
* WIP Multi-Big M * Deleting mbm copy.jl as it's not needed anymore * Added code coverage. * More test cases * Testing MBM type * working ver * Code coverage 100% for mbm.jl file. * Added comments, full code coverage on mbm.jl file TO BE TESTED * Code coverage tested to be 100%. Further documentation for mbm test_replace_variables_in_constraint() testcase hack. * Delete src/codecoverage.jl deleting prior to merge to main branch * Update datatypes.jl * Updated documentation, JuMP. calls. * More JuMP. additions * Works with AbstractVariableRef * Made _copy_variable function * _copy_variable function added. VariableRef changed to AbstractVariableRef. * Updated _copy_variable() * Update before deleting * Initial working solution. To be tested. * . * . * addition of option for user to specify final reformulation technique * working copy * test file added. * working to fix nested disjunctions * solve.jl works * 80 character limit * . * removed old code from constraints.jl * removed spacing edits. * . * . * working version with tests * Refactor variable handling and add model copying function Refactor variable creation and add copy_model_and_gdp_data function. * Add tests for copy_model_and_gdp_data function * copy_model integratred into cutting planes. keyword arguements added * added remapping functions for specific gdpdata field. created more tests accordingly. addition of test to check that model instances do not affect each other. * added solving test in modeljl for copy_model * Updated tests. * update tests * modified tests to check keyword arguments * updated cutting_planes function with keyword arguements * update readme * additional comments. unified spacing. * more spacing fixed * modified tests
1 parent 5709764 commit dc17a6e

File tree

7 files changed

+359
-4
lines changed

7 files changed

+359
-4
lines changed

‎README.md‎

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -184,6 +184,15 @@ The following reformulation methods are currently supported:
184184

185185
All variables must be included in exactly one partition. For manual partitioning, ensure each variable appears in exactly one group. For automatic partitioning, variables are divided as evenly as possible among the specified number of partitions.
186186

187+
6. [Cutting Planes](https://pubsonline.informs.org/doi/10.1287/ijoc.2015.0669): This method iteratively generates cutting planes using a separation problem and a relaxed Big-M formulation, then applies a final reformulation method. The `cutting_planes` struct is created with the following arguments:
188+
189+
- `optimizer`: Optimizer to use when solving the separation and relaxed Big-M subproblems. This is a required value.
190+
- `max_iter`: Maximum number of cutting plane iterations. Default: `3`.
191+
- `seperation_tolerance`: Convergence tolerance for the separation problem objective. Default: `1e-6`.
192+
- `final_reform_method`: Reformulation method to apply after cutting plane iterations. Default: `BigM()`.
193+
- `M_value`: Big-M value to use in the relaxed Big-M reformulation during iterations. Default: `1e9`.
194+
195+
187196

188197
## Release Notes
189198

‎src/DisjunctiveProgramming.jl‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -20,6 +20,7 @@ include("constraints.jl")
2020
include("macros.jl")
2121
include("reformulate.jl")
2222
include("bigm.jl")
23+
include("cuttingplanes.jl")
2324
include("hull.jl")
2425
include("mbm.jl")
2526
include("indicator.jl")

‎src/cuttingplanes.jl‎

Lines changed: 149 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,149 @@
1+
function reformulate_model(
2+
model::JuMP.AbstractModel,
3+
method::cutting_planes
4+
)
5+
_clear_reformulations(model)
6+
var_type = JuMP.variable_ref_type(model)
7+
obj = objective_function(model)
8+
sense = objective_sense(model)
9+
10+
#Creation of seperation (SEP) and relaxed big M model (rBM).
11+
SEP, sep_ref_map, _ = copy_gdp_model(model)
12+
rBM, rBM_ref_map, _ = copy_gdp_model(model)
13+
reformulate_model(rBM, BigM(method.M_value))
14+
reformulate_model(SEP, Hull())
15+
main_to_SEP_map = Dict(v => sep_ref_map[v] for v in all_variables(model))
16+
main_to_rBM_map = Dict(v => rBM_ref_map[v] for v in all_variables(model))
17+
JuMP.set_optimizer(SEP, method.optimizer)
18+
JuMP.set_optimizer(rBM, method.optimizer)
19+
JuMP.set_silent(rBM)
20+
JuMP.set_silent(SEP)
21+
JuMP.relax_integrality(rBM)
22+
JuMP.relax_integrality(SEP)
23+
JuMP.@objective(rBM, sense,
24+
_replace_variables_in_constraint(obj, main_to_rBM_map)
25+
)
26+
#Mapping of variables between models.
27+
rBM_to_SEP_map = Dict{var_type, var_type}()
28+
SEP_to_rBM_map = Dict{var_type, var_type}()
29+
for (var, rBM_var) in main_to_rBM_map
30+
SEP_var = main_to_SEP_map[var]
31+
rBM_to_SEP_map[rBM_var] = SEP_var
32+
SEP_to_rBM_map[SEP_var] = rBM_var
33+
end
34+
35+
#Main cutting planes loop.
36+
i = 1
37+
sep_obj = Inf
38+
while i <= method.max_iter && sep_obj > method.seperation_tolerance
39+
rBM_sol = _solve_rBM(rBM)
40+
SEP_sol = _solve_SEP(SEP, rBM, rBM_sol, SEP_to_rBM_map, rBM_to_SEP_map)
41+
sep_obj = objective_value(SEP)
42+
_cutting_planes(model, rBM, main_to_rBM_map,
43+
main_to_SEP_map, rBM_sol, SEP_sol
44+
)
45+
i += 1
46+
end
47+
48+
#Final reformulation with added cutting planes.
49+
reformulate_model(model, method.final_reform_method)
50+
return
51+
end
52+
53+
function _solve_rBM(
54+
rBM::M,
55+
) where {M <: JuMP.AbstractModel}
56+
T = JuMP.value_type(M)
57+
optimize!(rBM, ignore_optimize_hook = true)
58+
rBM_vars = JuMP.all_variables(rBM)
59+
60+
#Solution to be passed to SEP model.
61+
sol = Dict{JuMP.AbstractVariableRef,T}(var => zero(T) for var in rBM_vars)
62+
for rBM_var in rBM_vars
63+
sol[rBM_var] = JuMP.value(rBM_var)
64+
end
65+
return sol
66+
end
67+
68+
function _solve_SEP(
69+
SEP::M,
70+
rBM::M,
71+
rBM_sol::Dict{<:JuMP.AbstractVariableRef,T},
72+
SEP_to_rBM_map::Dict{<:JuMP.AbstractVariableRef,<:JuMP.AbstractVariableRef},
73+
rBM_to_SEP_map::Dict{<:JuMP.AbstractVariableRef,<:JuMP.AbstractVariableRef}
74+
) where {M <: JuMP.AbstractModel, T <: Number}
75+
76+
SEP_vars = [rBM_to_SEP_map[rBM_var] for rBM_var in JuMP.all_variables(rBM)]
77+
78+
#Modified objective function for SEP.
79+
obj_expr = sum(
80+
(SEP_var - rBM_sol[SEP_to_rBM_map[SEP_var]])^2 for SEP_var in SEP_vars
81+
)
82+
JuMP.@objective(SEP, Min, obj_expr)
83+
optimize!(SEP, ignore_optimize_hook = true)
84+
85+
#Solution to be used in cutting plane generation.
86+
sol = Dict{JuMP.AbstractVariableRef, T}(var => zero(T) for var in SEP_vars)
87+
for SEP_var in SEP_vars
88+
sol[SEP_var] = JuMP.value(SEP_var)
89+
end
90+
return sol
91+
end
92+
93+
function _cutting_planes(
94+
model::M,
95+
rBM::M,
96+
main_to_rBM_map::Dict{<:JuMP.AbstractVariableRef,<:JuMP.AbstractVariableRef},
97+
main_to_SEP_map::Dict{<:JuMP.AbstractVariableRef,<:JuMP.AbstractVariableRef},
98+
rBM_sol::Dict{<:JuMP.AbstractVariableRef,T},
99+
SEP_sol::Dict{<:JuMP.AbstractVariableRef,T},
100+
) where {M <: JuMP.AbstractModel, T <: Number}
101+
main_vars = JuMP.all_variables(model)
102+
103+
#Cutting plane generation
104+
ξ_sep = Dict{JuMP.AbstractVariableRef,T}(var =>zero(T) for var in main_vars)
105+
for var in main_vars
106+
ξ_sep[var] = 2*(SEP_sol[main_to_SEP_map[var]]
107+
-rBM_sol[main_to_rBM_map[var]]
108+
)
109+
end
110+
#Cutting plane added to main model.
111+
main_cut = JuMP.@expression(model,
112+
sum(ξ_sep[var]*(var - SEP_sol[main_to_SEP_map[var]])
113+
for var in main_vars
114+
)
115+
)
116+
#Cutting plane added to rBM
117+
rBM_cut = _replace_variables_in_constraint(main_cut, main_to_rBM_map)
118+
JuMP.@constraint(model, main_cut >= 0.0)
119+
JuMP.@constraint(rBM, rBM_cut >= 0.0)
120+
end
121+
122+
################################################################################
123+
# ERROR MESSAGES
124+
################################################################################
125+
126+
function reformulate_model(::M, ::cutting_planes) where {M}
127+
error("reformulate_model not implemented for model type `$(M)`.")
128+
end
129+
130+
function _solve_rBM(::M) where {M}
131+
error("_solve_rBM not implemented for model type `$(M)`.")
132+
end
133+
134+
function _solve_SEP(::M, ::N, ::H, ::S, ::R) where {M, N, H, S, R}
135+
error("_solve_SEP not implemented for argument types:\n
136+
SEP: `$(M)`, rBM: `$(N)`,\n
137+
rBM_sol: `$(H)`,\n
138+
SEP_to_rBM_map: `$(S)`,\n
139+
rBM_to_SEP_map: `$(R)`.")
140+
end
141+
142+
function _cutting_planes(::M, ::N, ::H, ::S, ::R, ::T) where {M, N, H, S, R, T}
143+
error("_cutting_planes not implemented for argument types: \n
144+
model: `$(M)`, rBM: `$(N)`,\n
145+
main_to_rBM_map: `$(H)`, main_to_SEP_map:
146+
`$(S)`,\n
147+
rBM_sol: `$(R)`,\n
148+
SEP_sol: `$(T)`.")
149+
end

‎src/datatypes.jl‎

Lines changed: 30 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -431,6 +431,36 @@ mutable struct _Hull{V <: JuMP.AbstractVariableRef, T} <: AbstractReformulationM
431431
end
432432
end
433433

434+
"""
435+
cutting_planes{O} <: AbstractReformulationMethod
436+
437+
A type for using the cutting planes approach for disjunctive constraints.
438+
439+
**Fields**
440+
- `optimizer::O`: Optimizer to use when solving mini-models (required).
441+
- `max_iter::Int`: Number of iterations (default = `3`).
442+
- `seperation_tolerance::Float64`: Tolerance for the separation problem (default = `1e-6`).
443+
- `final_reform_method::AbstractReformulationMethod`: Final reformulation
444+
method to use after cutting planes (default = `BigM()`).
445+
- `M_value::Float64`: Big-M value to use in the final reformulation (default = `1e9`).
446+
"""
447+
struct cutting_planes{O} <: AbstractReformulationMethod
448+
optimizer::O;
449+
max_iter::Int
450+
seperation_tolerance::Float64
451+
final_reform_method::AbstractReformulationMethod
452+
M_value::Float64
453+
function cutting_planes(
454+
optimizer::O;
455+
max_iter::Int = 3,
456+
seperation_tolerance::Float64 = 1e-6,
457+
final_reform_method = BigM(),
458+
M_value::Float64 = 1e9
459+
) where {O}
460+
new{O}(optimizer, max_iter, seperation_tolerance, final_reform_method, M_value)
461+
end
462+
end
463+
434464
"""
435465
PSplit <: AbstractReformulationMethod
436466

‎test/constraints/cuttingplanes.jl‎

Lines changed: 155 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,155 @@
1+
using HiGHS
2+
3+
function test_cutting_planes_datatype()
4+
method = cutting_planes(HiGHS.Optimizer)
5+
@test method.optimizer == HiGHS.Optimizer
6+
@test method.max_iter == 3
7+
@test method.seperation_tolerance == 1e-6
8+
@test method.final_reform_method isa BigM
9+
@test method.M_value == 1e9
10+
11+
method = cutting_planes(HiGHS.Optimizer;max_iter=10,
12+
seperation_tolerance=1e-4, final_reform_method=Indicator(), M_value=1e6
13+
)
14+
@test method.max_iter == 10
15+
@test method.seperation_tolerance == 1e-4
16+
@test method.final_reform_method isa Indicator
17+
@test method.M_value == 1e6
18+
end
19+
20+
function test_solve_rBM()
21+
rBM = JuMP.Model(HiGHS.Optimizer)
22+
@variable(rBM, 0 <= x <= 100)
23+
@variable(rBM, 0 <= y[1:2] <= 1)
24+
@constraint(rBM, x <= 3 + 100(1 - y[1]))
25+
@constraint(rBM, x <= 4 + 100(1 - y[2]))
26+
@constraint(rBM, y[1] + y[2] == 1)
27+
@objective(rBM, Max, x)
28+
29+
solutions = DP._solve_rBM(rBM)
30+
@test solutions[x] == 53.5
31+
@test solutions[y[1]] == 0.495
32+
@test solutions[y[2]] == 0.505
33+
34+
@test_throws ErrorException DP._solve_rBM(Dict())
35+
end
36+
37+
function test_solve_SEP()
38+
model = GDPModel()
39+
@variable(model, 0 <= x <= 100)
40+
@variable(model, Y[1:2], Logical)
41+
@constraint(model, x <= 3, Disjunct(Y[1]))
42+
@constraint(model, x <= 4, Disjunct(Y[2]))
43+
@disjunction(model, [Y[1], Y[2]])
44+
@objective(model, Max, x)
45+
var_type = JuMP.variable_ref_type(model)
46+
method = cutting_planes(HiGHS.Optimizer)
47+
obj = objective_function(model)
48+
sense = objective_sense(model)
49+
SEP, sep_ref_map, _ = DP.copy_gdp_model(model)
50+
rBM, rBM_ref_map, _ = DP.copy_gdp_model(model)
51+
DP.reformulate_model(rBM, DP.BigM(method.M_value))
52+
DP.reformulate_model(SEP, DP.Hull())
53+
main_to_SEP_map = Dict(v => sep_ref_map[v] for v in all_variables(model))
54+
main_to_rBM_map = Dict(v => rBM_ref_map[v] for v in all_variables(model))
55+
JuMP.set_optimizer(SEP, method.optimizer)
56+
JuMP.set_optimizer(rBM, method.optimizer)
57+
JuMP.set_silent(rBM)
58+
JuMP.set_silent(SEP)
59+
JuMP.relax_integrality(rBM)
60+
JuMP.relax_integrality(SEP)
61+
JuMP.@objective(rBM, sense,
62+
DP._replace_variables_in_constraint(obj, main_to_rBM_map)
63+
)
64+
rBM_to_SEP_map = Dict{var_type, var_type}()
65+
SEP_to_rBM_map = Dict{var_type, var_type}()
66+
for (var, rBM_var) in main_to_rBM_map
67+
SEP_var = main_to_SEP_map[var]
68+
rBM_to_SEP_map[rBM_var] = SEP_var
69+
SEP_to_rBM_map[SEP_var] = rBM_var
70+
end
71+
rBM_sol = DP._solve_rBM(rBM)
72+
SEP_sol = DP._solve_SEP(SEP, rBM, rBM_sol, SEP_to_rBM_map, rBM_to_SEP_map)
73+
@test length(SEP_sol) == length(rBM_sol)
74+
@test SEP_sol[rBM_to_SEP_map[main_to_rBM_map[x]]] ≈ 4.0
75+
76+
@test_throws ErrorException DP._solve_SEP(SEP, rBM, rBM_sol, SEP_to_rBM_map
77+
, "not a dict"
78+
)
79+
end
80+
81+
function test_cutting_planes()
82+
model = GDPModel()
83+
@variable(model, 0 <= x <= 100)
84+
@variable(model, Y[1:2], Logical)
85+
@constraint(model, x <= 3, Disjunct(Y[1]))
86+
@constraint(model, x <= 4, Disjunct(Y[2]))
87+
@disjunction(model, [Y[1], Y[2]])
88+
@objective(model, Max, x)
89+
var_type = JuMP.variable_ref_type(model)
90+
method = cutting_planes(HiGHS.Optimizer)
91+
obj = objective_function(model)
92+
sense = objective_sense(model)
93+
SEP, sep_ref_map, _ = DP.copy_gdp_model(model)
94+
rBM, rBM_ref_map, _ = DP.copy_gdp_model(model)
95+
DP.reformulate_model(rBM, DP.BigM(method.M_value))
96+
DP.reformulate_model(SEP, DP.Hull())
97+
main_to_SEP_map = Dict(v => sep_ref_map[v] for v in all_variables(model))
98+
main_to_rBM_map = Dict(v => rBM_ref_map[v] for v in all_variables(model))
99+
JuMP.set_optimizer(SEP, method.optimizer)
100+
JuMP.set_optimizer(rBM, method.optimizer)
101+
JuMP.set_silent(rBM)
102+
JuMP.set_silent(SEP)
103+
JuMP.relax_integrality(rBM)
104+
JuMP.relax_integrality(SEP)
105+
JuMP.@objective(rBM, sense,
106+
DP._replace_variables_in_constraint(obj, main_to_rBM_map)
107+
)
108+
rBM_to_SEP_map = Dict{var_type, var_type}()
109+
SEP_to_rBM_map = Dict{var_type, var_type}()
110+
for (var, rBM_var) in main_to_rBM_map
111+
SEP_var = main_to_SEP_map[var]
112+
rBM_to_SEP_map[rBM_var] = SEP_var
113+
SEP_to_rBM_map[SEP_var] = rBM_var
114+
end
115+
rBM_sol = DP._solve_rBM(rBM)
116+
SEP_sol = DP._solve_SEP(SEP, rBM, rBM_sol, SEP_to_rBM_map, rBM_to_SEP_map)
117+
DP._cutting_planes(model, rBM, main_to_rBM_map, main_to_SEP_map, rBM_sol, SEP_sol)
118+
119+
rBM_sol = DP._solve_rBM(rBM)
120+
SEP_sol = DP._solve_SEP(SEP, rBM, rBM_sol, SEP_to_rBM_map, rBM_to_SEP_map)
121+
122+
@test rBM_sol[main_to_rBM_map[x]] ≈ 4.0
123+
@test SEP_sol[rBM_to_SEP_map[main_to_rBM_map[x]]] ≈ 4.0 atol=1e-3
124+
125+
@test_throws ErrorException DP._cutting_planes(model, rBM, main_to_rBM_map,
126+
main_to_SEP_map, rBM_sol, "not a dict"
127+
)
128+
end
129+
130+
function test_reformulate_model()
131+
model = GDPModel()
132+
@variable(model, 0 <= x[1:4] <= 100)
133+
@variable(model, Y[1:2], Logical)
134+
@constraint(model, x[1] + x[2] <= 3, Disjunct(Y[1]))
135+
@constraint(model, x[3] + x[4] <= 4, Disjunct(Y[2]))
136+
@disjunction(model, [Y[1], Y[2]])
137+
@objective(model, Max, x[1] + x[2])
138+
139+
method = cutting_planes(HiGHS.Optimizer)
140+
DP.reformulate_model(model, method)
141+
num_con = length(
142+
JuMP.all_constraints(model; include_variable_in_set_constraints = false)
143+
)
144+
@test num_con == 4
145+
@test_throws ErrorException DP.reformulate_model(42, method)
146+
end
147+
148+
149+
@testset "Cutting Planes" begin
150+
test_cutting_planes_datatype()
151+
test_solve_rBM()
152+
test_solve_SEP()
153+
test_cutting_planes()
154+
test_reformulate_model()
155+
end

‎test/runtests.jl‎

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -17,8 +17,9 @@ include("constraints/indicator.jl")
1717
include("constraints/mbm.jl")
1818
include("constraints/bigm.jl")
1919
include("constraints/psplit.jl")
20+
include("constraints/cuttingplanes.jl")
2021
include("constraints/hull.jl")
2122
include("constraints/fallback.jl")
2223
include("constraints/disjunction.jl")
2324
include("print.jl")
24-
include("solve.jl")
25+
include("solve.jl")

0 commit comments

Comments
 (0)