1# ___________________________________________________________________________ 2# 3# Pyomo: Python Optimization Modeling Objects 4# Copyright 2017 National Technology and Engineering Solutions of Sandia, LLC 5# Under the terms of Contract DE-NA0003525 with National Technology and 6# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain 7# rights in this software. 8# This software is distributed under the 3-clause BSD License. 9# ___________________________________________________________________________ 10 11"""Tests for the GDPopt solver plugin.""" 12import logging 13from math import fabs 14from os.path import join, normpath 15 16from io import StringIO 17 18import pyomo.common.unittest as unittest 19from pyomo.common.log import LoggingIntercept 20from pyomo.common.collections import Bunch 21from pyomo.common.fileutils import import_file 22from pyomo.contrib.gdpopt.GDPopt import GDPoptSolver 23from pyomo.contrib.gdpopt.data_class import GDPoptSolveData 24from pyomo.contrib.gdpopt.mip_solve import solve_linear_GDP 25from pyomo.contrib.gdpopt.util import is_feasible, time_code 26from pyomo.environ import ConcreteModel, Objective, SolverFactory, Var, value, Integers, Block, Constraint, maximize 27from pyomo.gdp import Disjunct, Disjunction 28from pyomo.contrib.mcpp.pyomo_mcpp import mcpp_available 29from pyomo.opt import TerminationCondition 30 31from pyomo.common.fileutils import PYOMO_ROOT_DIR 32exdir = normpath(join(PYOMO_ROOT_DIR, 'examples', 'gdp')) 33 34mip_solver = 'glpk' 35nlp_solver = 'ipopt' 36global_nlp_solver = 'baron' 37global_nlp_solver_args = dict() 38minlp_solver = 'baron' 39LOA_solvers = (mip_solver, nlp_solver) 40GLOA_solvers = (mip_solver, global_nlp_solver, minlp_solver) 41LOA_solvers_available = all(SolverFactory(s).available() for s in LOA_solvers) 42GLOA_solvers_available = all(SolverFactory(s).available() for s in GLOA_solvers) 43license_available = SolverFactory(global_nlp_solver).license_is_valid() if GLOA_solvers_available else False 44 45 46class TestGDPoptUnit(unittest.TestCase): 47 """Real unit tests for GDPopt""" 48 49 @unittest.skipUnless(SolverFactory(mip_solver).available(), "MIP solver not available") 50 def test_solve_linear_GDP_unbounded(self): 51 m = ConcreteModel() 52 m.GDPopt_utils = Block() 53 m.x = Var(bounds=(-1, 10)) 54 m.y = Var(bounds=(2, 3)) 55 m.z = Var() 56 m.d = Disjunction(expr=[ 57 [m.x + m.y >= 5], [m.x - m.y <= 3] 58 ]) 59 m.o = Objective(expr=m.z) 60 m.GDPopt_utils.variable_list = [m.x, m.y, m.z] 61 m.GDPopt_utils.disjunct_list = [m.d._autodisjuncts[0], m.d._autodisjuncts[1]] 62 output = StringIO() 63 with LoggingIntercept(output, 'pyomo.contrib.gdpopt', logging.WARNING): 64 solver_data = GDPoptSolveData() 65 solver_data.timing = Bunch() 66 with time_code(solver_data.timing, 'main', is_main_timer=True): 67 solve_linear_GDP(m, solver_data, GDPoptSolver.CONFIG(dict(mip_solver=mip_solver, strategy='LOA'))) 68 self.assertIn("Linear GDP was unbounded. Resolving with arbitrary bound values", 69 output.getvalue().strip()) 70 71 @unittest.skipUnless(SolverFactory(mip_solver).available(), "MIP solver not available") 72 def test_solve_lp(self): 73 m = ConcreteModel() 74 m.x = Var(bounds=(-5, 5)) 75 m.c = Constraint(expr=m.x >= 1) 76 m.o = Objective(expr=m.x) 77 output = StringIO() 78 with LoggingIntercept(output, 'pyomo.contrib.gdpopt', logging.INFO): 79 SolverFactory('gdpopt').solve(m, mip_solver=mip_solver, strategy='LOA') 80 self.assertIn("Your model is an LP (linear program).", 81 output.getvalue().strip()) 82 self.assertAlmostEqual(value(m.o.expr), 1) 83 84 @unittest.skipUnless(SolverFactory(nlp_solver).available(), 'NLP solver not available') 85 def test_solve_nlp(self): 86 m = ConcreteModel() 87 m.x = Var(bounds=(-5, 5)) 88 m.c = Constraint(expr=m.x >= 1) 89 m.o = Objective(expr=m.x ** 2) 90 output = StringIO() 91 with LoggingIntercept(output, 'pyomo.contrib.gdpopt', logging.INFO): 92 SolverFactory('gdpopt').solve(m, nlp_solver=nlp_solver, strategy='LOA') 93 self.assertIn("Your model is an NLP (nonlinear program).", 94 output.getvalue().strip()) 95 self.assertAlmostEqual(value(m.o.expr), 1) 96 97 @unittest.skipUnless(SolverFactory(mip_solver).available(), "MIP solver not available") 98 def test_solve_constant_obj(self): 99 m = ConcreteModel() 100 m.x = Var(bounds=(-5, 5)) 101 m.c = Constraint(expr=m.x >= 1) 102 m.o = Objective(expr=1) 103 output = StringIO() 104 with LoggingIntercept(output, 'pyomo.contrib.gdpopt', logging.INFO): 105 SolverFactory('gdpopt').solve(m, mip_solver=mip_solver, strategy='LOA') 106 self.assertIn("Your model is an LP (linear program).", 107 output.getvalue().strip()) 108 self.assertAlmostEqual(value(m.o.expr), 1) 109 110 @unittest.skipUnless(SolverFactory(nlp_solver).available(), 'NLP solver not available') 111 def test_no_objective(self): 112 m = ConcreteModel() 113 m.x = Var(bounds=(-5, 5)) 114 m.c = Constraint(expr=m.x ** 2 >= 1) 115 output = StringIO() 116 with LoggingIntercept(output, 'pyomo.contrib.gdpopt', logging.WARNING): 117 SolverFactory('gdpopt').solve(m, nlp_solver=nlp_solver, strategy='LOA') 118 self.assertIn("Model has no active objectives. Adding dummy objective.", 119 output.getvalue().strip()) 120 121 def test_multiple_objectives(self): 122 m = ConcreteModel() 123 m.x = Var() 124 m.o = Objective(expr=m.x) 125 m.o2 = Objective(expr=m.x + 1) 126 with self.assertRaisesRegex(ValueError, "Model has multiple active objectives"): 127 SolverFactory('gdpopt').solve(m, strategy='LOA') 128 129 def test_is_feasible_function(self): 130 m = ConcreteModel() 131 m.x = Var(bounds=(0, 3), initialize=2) 132 m.c = Constraint(expr=m.x == 2) 133 self.assertTrue(is_feasible(m, GDPoptSolver.CONFIG(dict(strategy='LOA')))) 134 135 m.c2 = Constraint(expr=m.x <= 1) 136 self.assertFalse(is_feasible(m, GDPoptSolver.CONFIG(dict(strategy='LOA')))) 137 138 m = ConcreteModel() 139 m.x = Var(bounds=(0, 3), initialize=2) 140 m.c = Constraint(expr=m.x >= 5) 141 self.assertFalse(is_feasible(m, GDPoptSolver.CONFIG(dict(strategy='LOA')))) 142 143 m = ConcreteModel() 144 m.x = Var(bounds=(3, 3), initialize=2) 145 self.assertFalse(is_feasible(m, GDPoptSolver.CONFIG(dict(strategy='LOA')))) 146 147 m = ConcreteModel() 148 m.x = Var(bounds=(0, 1), initialize=2) 149 self.assertFalse(is_feasible(m, GDPoptSolver.CONFIG(dict(strategy='LOA')))) 150 151 m = ConcreteModel() 152 m.x = Var(bounds=(0, 1), initialize=2) 153 m.d = Disjunct() 154 with self.assertRaisesRegex(NotImplementedError, "Found active disjunct"): 155 is_feasible(m, GDPoptSolver.CONFIG(dict(strategy='LOA'))) 156 157 158@unittest.skipIf(not LOA_solvers_available, 159 "Required subsolvers %s are not available" 160 % (LOA_solvers,)) 161class TestGDPopt(unittest.TestCase): 162 """Tests for the GDPopt solver plugin.""" 163 164 def test_infeasible_GDP(self): 165 """Test for infeasible GDP.""" 166 m = ConcreteModel() 167 m.x = Var(bounds=(0, 2)) 168 m.d = Disjunction(expr=[ 169 [m.x ** 2 >= 3, m.x >= 3], 170 [m.x ** 2 <= -1, m.x <= -1]]) 171 m.o = Objective(expr=m.x) 172 output = StringIO() 173 with LoggingIntercept(output, 'pyomo.contrib.gdpopt', logging.WARNING): 174 SolverFactory('gdpopt').solve( 175 m, strategy='LOA', 176 mip_solver=mip_solver, 177 nlp_solver=nlp_solver) 178 self.assertIn("Set covering problem was infeasible.", 179 output.getvalue().strip()) 180 181 def test_GDP_nonlinear_objective(self): 182 m = ConcreteModel() 183 m.x = Var(bounds=(-1, 10)) 184 m.y = Var(bounds=(2, 3)) 185 m.d = Disjunction(expr=[ 186 [m.x + m.y >= 5], [m.x - m.y <= 3] 187 ]) 188 m.o = Objective(expr=m.x ** 2) 189 SolverFactory('gdpopt').solve( 190 m, strategy='LOA', 191 mip_solver=mip_solver, 192 nlp_solver=nlp_solver 193 ) 194 self.assertAlmostEqual(value(m.o), 0) 195 196 m = ConcreteModel() 197 m.x = Var(bounds=(-1, 10)) 198 m.y = Var(bounds=(2, 3)) 199 m.d = Disjunction(expr=[ 200 [m.x + m.y >= 5], [m.x - m.y <= 3] 201 ]) 202 m.o = Objective(expr=-m.x ** 2, sense=maximize) 203 SolverFactory('gdpopt').solve( 204 m, strategy='LOA', 205 mip_solver=mip_solver, 206 nlp_solver=nlp_solver 207 ) 208 self.assertAlmostEqual(value(m.o), 0) 209 210 def test_LOA_8PP_default_init(self): 211 """Test logic-based outer approximation with 8PP.""" 212 exfile = import_file( 213 join(exdir, 'eight_process', 'eight_proc_model.py')) 214 eight_process = exfile.build_eight_process_flowsheet() 215 SolverFactory('gdpopt').solve( 216 eight_process, strategy='LOA', 217 mip_solver=mip_solver, 218 nlp_solver=nlp_solver, 219 tee=False) 220 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 221 222 @unittest.skipUnless(SolverFactory('gams').available(exception_flag=False), 'GAMS solver not available') 223 def test_LOA_8PP_gams_solver(self): 224 # Make sure that the duals are still correct 225 exfile = import_file( 226 join(exdir, 'eight_process', 'eight_proc_model.py')) 227 eight_process = exfile.build_eight_process_flowsheet() 228 SolverFactory('gdpopt').solve( 229 eight_process, strategy='LOA', 230 mip_solver=mip_solver, 231 nlp_solver='gams', 232 max_slack=0, 233 tee=False) 234 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 235 236 def test_LOA_8PP_force_NLP(self): 237 exfile = import_file( 238 join(exdir, 'eight_process', 'eight_proc_model.py')) 239 eight_process = exfile.build_eight_process_flowsheet() 240 SolverFactory('gdpopt').solve( 241 eight_process, strategy='LOA', 242 mip_solver=mip_solver, 243 nlp_solver=nlp_solver, 244 force_subproblem_nlp=True, 245 tee=False) 246 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 247 248 def test_LOA_strip_pack_default_init(self): 249 """Test logic-based outer approximation with strip packing.""" 250 exfile = import_file( 251 join(exdir, 'strip_packing', 'strip_packing_concrete.py')) 252 strip_pack = exfile.build_rect_strip_packing_model() 253 SolverFactory('gdpopt').solve( 254 strip_pack, strategy='LOA', 255 mip_solver=mip_solver, 256 nlp_solver=nlp_solver) 257 self.assertTrue( 258 fabs(value(strip_pack.total_length.expr) - 11) <= 1E-2) 259 260 @unittest.category('expensive') 261 def test_LOA_constrained_layout_default_init(self): 262 """Test LOA with constrained layout.""" 263 exfile = import_file( 264 join(exdir, 'constrained_layout', 'cons_layout_model.py')) 265 cons_layout = exfile.build_constrained_layout_model() 266 SolverFactory('gdpopt').solve( 267 cons_layout, strategy='LOA', 268 mip_solver=mip_solver, 269 nlp_solver=nlp_solver, 270 iterlim=120, 271 max_slack=5, # problem is convex, so can decrease slack 272 ) 273 objective_value = value(cons_layout.min_dist_cost.expr) 274 self.assertTrue( 275 fabs(objective_value - 41573) <= 200, 276 "Objective value of %s instead of 41573" % objective_value) 277 278 def test_LOA_8PP_maxBinary(self): 279 """Test logic-based OA with max_binary initialization.""" 280 exfile = import_file( 281 join(exdir, 'eight_process', 'eight_proc_model.py')) 282 eight_process = exfile.build_eight_process_flowsheet() 283 SolverFactory('gdpopt').solve( 284 eight_process, strategy='LOA', init_strategy='max_binary', 285 mip_solver=mip_solver, 286 nlp_solver=nlp_solver) 287 288 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 289 290 def test_LOA_strip_pack_maxBinary(self): 291 """Test LOA with strip packing using max_binary initialization.""" 292 exfile = import_file( 293 join(exdir, 'strip_packing', 'strip_packing_concrete.py')) 294 strip_pack = exfile.build_rect_strip_packing_model() 295 SolverFactory('gdpopt').solve( 296 strip_pack, strategy='LOA', init_strategy='max_binary', 297 mip_solver=mip_solver, 298 nlp_solver=nlp_solver) 299 self.assertTrue( 300 fabs(value(strip_pack.total_length.expr) - 11) <= 1E-2) 301 302 def test_LOA_8PP_fixed_disjuncts(self): 303 """Test LOA with 8PP using fixed disjuncts initialization.""" 304 exfile = import_file( 305 join(exdir, 'eight_process', 'eight_proc_model.py')) 306 eight_process = exfile.build_eight_process_flowsheet() 307 initialize = [ 308 # Use units 1, 4, 7, 8 309 eight_process.use_unit_1or2.disjuncts[0], 310 eight_process.use_unit_3ornot.disjuncts[1], 311 eight_process.use_unit_4or5ornot.disjuncts[0], 312 eight_process.use_unit_6or7ornot.disjuncts[1], 313 eight_process.use_unit_8ornot.disjuncts[0] 314 ] 315 for disj in eight_process.component_data_objects(Disjunct): 316 if disj in initialize: 317 disj.indicator_var.set_value(1) 318 else: 319 disj.indicator_var.set_value(0) 320 SolverFactory('gdpopt').solve( 321 eight_process, strategy='LOA', init_strategy='fix_disjuncts', 322 mip_solver=mip_solver, 323 nlp_solver=nlp_solver) 324 325 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 326 327 def test_LOA_custom_disjuncts(self): 328 """Test logic-based OA with custom disjuncts initialization.""" 329 exfile = import_file( 330 join(exdir, 'eight_process', 'eight_proc_model.py')) 331 eight_process = exfile.build_eight_process_flowsheet() 332 initialize = [ 333 # Use units 1, 4, 7, 8 334 [eight_process.use_unit_1or2.disjuncts[0], 335 eight_process.use_unit_3ornot.disjuncts[1], 336 eight_process.use_unit_4or5ornot.disjuncts[0], 337 eight_process.use_unit_6or7ornot.disjuncts[1], 338 eight_process.use_unit_8ornot.disjuncts[0]], 339 # Use units 2, 4, 6, 8 340 [eight_process.use_unit_1or2.disjuncts[1], 341 eight_process.use_unit_3ornot.disjuncts[1], 342 eight_process.use_unit_4or5ornot.disjuncts[0], 343 eight_process.use_unit_6or7ornot.disjuncts[0], 344 eight_process.use_unit_8ornot.disjuncts[0]] 345 ] 346 347 def assert_correct_disjuncts_active(nlp_model, solve_data): 348 if solve_data.master_iteration >= 1: 349 return # only checking initialization 350 iter_num = solve_data.nlp_iteration 351 disjs_should_be_active = initialize[iter_num - 1] 352 for orig_disj, soln_disj in zip( 353 solve_data.original_model.GDPopt_utils.disjunct_list, 354 nlp_model.GDPopt_utils.disjunct_list 355 ): 356 if orig_disj in disjs_should_be_active: 357 self.assertTrue(soln_disj.indicator_var.value == 1) 358 359 SolverFactory('gdpopt').solve( 360 eight_process, strategy='LOA', init_strategy='custom_disjuncts', 361 custom_init_disjuncts=initialize, 362 mip_solver=mip_solver, 363 nlp_solver=nlp_solver, 364 call_after_subproblem_feasible=assert_correct_disjuncts_active) 365 366 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 367 368 369@unittest.skipIf(not LOA_solvers_available, 370 "Required subsolvers %s are not available" 371 % (LOA_solvers,)) 372class TestGDPoptRIC(unittest.TestCase): 373 """Tests for the GDPopt solver plugin.""" 374 375 def test_infeasible_GDP(self): 376 """Test for infeasible GDP.""" 377 m = ConcreteModel() 378 m.x = Var(bounds=(0, 2)) 379 m.d = Disjunction(expr=[ 380 [m.x ** 2 >= 3, m.x >= 3], 381 [m.x ** 2 <= -1, m.x <= -1]]) 382 m.o = Objective(expr=m.x) 383 output = StringIO() 384 with LoggingIntercept(output, 'pyomo.contrib.gdpopt', logging.WARNING): 385 SolverFactory('gdpopt').solve( 386 m, strategy='RIC', 387 mip_solver=mip_solver, 388 nlp_solver=nlp_solver) 389 self.assertIn("Set covering problem was infeasible.", 390 output.getvalue().strip()) 391 392 def test_GDP_nonlinear_objective(self): 393 m = ConcreteModel() 394 m.x = Var(bounds=(-1, 10)) 395 m.y = Var(bounds=(2, 3)) 396 m.d = Disjunction(expr=[ 397 [m.x + m.y >= 5], [m.x - m.y <= 3] 398 ]) 399 m.o = Objective(expr=m.x ** 2) 400 SolverFactory('gdpopt').solve( 401 m, strategy='RIC', 402 mip_solver=mip_solver, 403 nlp_solver=nlp_solver 404 ) 405 self.assertAlmostEqual(value(m.o), 0) 406 407 m = ConcreteModel() 408 m.x = Var(bounds=(-1, 10)) 409 m.y = Var(bounds=(2, 3)) 410 m.d = Disjunction(expr=[ 411 [m.x + m.y >= 5], [m.x - m.y <= 3] 412 ]) 413 m.o = Objective(expr=-m.x ** 2, sense=maximize) 414 SolverFactory('gdpopt').solve( 415 m, strategy='RIC', 416 mip_solver=mip_solver, 417 nlp_solver=nlp_solver 418 ) 419 self.assertAlmostEqual(value(m.o), 0) 420 421 def test_RIC_8PP_default_init(self): 422 """Test logic-based outer approximation with 8PP.""" 423 exfile = import_file( 424 join(exdir, 'eight_process', 'eight_proc_model.py')) 425 eight_process = exfile.build_eight_process_flowsheet() 426 SolverFactory('gdpopt').solve( 427 eight_process, strategy='RIC', 428 mip_solver=mip_solver, 429 nlp_solver=nlp_solver, 430 tee=False) 431 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 432 433 @unittest.skipUnless(SolverFactory('gams').available(exception_flag=False), 'GAMS solver not available') 434 def test_RIC_8PP_gams_solver(self): 435 # Make sure that the duals are still correct 436 exfile = import_file( 437 join(exdir, 'eight_process', 'eight_proc_model.py')) 438 eight_process = exfile.build_eight_process_flowsheet() 439 SolverFactory('gdpopt').solve( 440 eight_process, strategy='RIC', 441 mip_solver=mip_solver, 442 nlp_solver='gams', 443 max_slack=0, 444 tee=False) 445 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 446 447 def test_RIC_8PP_force_NLP(self): 448 exfile = import_file( 449 join(exdir, 'eight_process', 'eight_proc_model.py')) 450 eight_process = exfile.build_eight_process_flowsheet() 451 SolverFactory('gdpopt').solve( 452 eight_process, strategy='RIC', 453 mip_solver=mip_solver, 454 nlp_solver=nlp_solver, 455 force_subproblem_nlp=True, 456 tee=False) 457 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 458 459 def test_RIC_strip_pack_default_init(self): 460 """Test logic-based outer approximation with strip packing.""" 461 exfile = import_file( 462 join(exdir, 'strip_packing', 'strip_packing_concrete.py')) 463 strip_pack = exfile.build_rect_strip_packing_model() 464 SolverFactory('gdpopt').solve( 465 strip_pack, strategy='RIC', 466 mip_solver=mip_solver, 467 nlp_solver=nlp_solver) 468 self.assertTrue( 469 fabs(value(strip_pack.total_length.expr) - 11) <= 1E-2) 470 471 @unittest.category('expensive') 472 def test_RIC_constrained_layout_default_init(self): 473 """Test RIC with constrained layout.""" 474 exfile = import_file( 475 join(exdir, 'constrained_layout', 'cons_layout_model.py')) 476 cons_layout = exfile.build_constrained_layout_model() 477 SolverFactory('gdpopt').solve( 478 cons_layout, strategy='RIC', 479 mip_solver=mip_solver, 480 nlp_solver=nlp_solver, 481 iterlim=120, 482 max_slack=5, # problem is convex, so can decrease slack 483 ) 484 objective_value = value(cons_layout.min_dist_cost.expr) 485 self.assertTrue( 486 fabs(objective_value - 41573) <= 200, 487 "Objective value of %s instead of 41573" % objective_value) 488 489 def test_RIC_8PP_maxBinary(self): 490 """Test logic-based OA with max_binary initialization.""" 491 exfile = import_file( 492 join(exdir, 'eight_process', 'eight_proc_model.py')) 493 eight_process = exfile.build_eight_process_flowsheet() 494 SolverFactory('gdpopt').solve( 495 eight_process, strategy='RIC', init_strategy='max_binary', 496 mip_solver=mip_solver, 497 nlp_solver=nlp_solver) 498 499 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 500 501 def test_RIC_strip_pack_maxBinary(self): 502 """Test RIC with strip packing using max_binary initialization.""" 503 exfile = import_file( 504 join(exdir, 'strip_packing', 'strip_packing_concrete.py')) 505 strip_pack = exfile.build_rect_strip_packing_model() 506 SolverFactory('gdpopt').solve( 507 strip_pack, strategy='RIC', init_strategy='max_binary', 508 mip_solver=mip_solver, 509 nlp_solver=nlp_solver) 510 self.assertTrue( 511 fabs(value(strip_pack.total_length.expr) - 11) <= 1E-2) 512 513 def test_RIC_8PP_fixed_disjuncts(self): 514 """Test RIC with 8PP using fixed disjuncts initialization.""" 515 exfile = import_file( 516 join(exdir, 'eight_process', 'eight_proc_model.py')) 517 eight_process = exfile.build_eight_process_flowsheet() 518 initialize = [ 519 # Use units 1, 4, 7, 8 520 eight_process.use_unit_1or2.disjuncts[0], 521 eight_process.use_unit_3ornot.disjuncts[1], 522 eight_process.use_unit_4or5ornot.disjuncts[0], 523 eight_process.use_unit_6or7ornot.disjuncts[1], 524 eight_process.use_unit_8ornot.disjuncts[0] 525 ] 526 for disj in eight_process.component_data_objects(Disjunct): 527 if disj in initialize: 528 disj.indicator_var.set_value(1) 529 else: 530 disj.indicator_var.set_value(0) 531 SolverFactory('gdpopt').solve( 532 eight_process, strategy='RIC', init_strategy='fix_disjuncts', 533 mip_solver=mip_solver, 534 nlp_solver=nlp_solver) 535 536 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 537 538 def test_RIC_custom_disjuncts(self): 539 """Test logic-based OA with custom disjuncts initialization.""" 540 exfile = import_file( 541 join(exdir, 'eight_process', 'eight_proc_model.py')) 542 eight_process = exfile.build_eight_process_flowsheet() 543 initialize = [ 544 # Use units 1, 4, 7, 8 545 [eight_process.use_unit_1or2.disjuncts[0], 546 eight_process.use_unit_3ornot.disjuncts[1], 547 eight_process.use_unit_4or5ornot.disjuncts[0], 548 eight_process.use_unit_6or7ornot.disjuncts[1], 549 eight_process.use_unit_8ornot.disjuncts[0]], 550 # Use units 2, 4, 6, 8 551 [eight_process.use_unit_1or2.disjuncts[1], 552 eight_process.use_unit_3ornot.disjuncts[1], 553 eight_process.use_unit_4or5ornot.disjuncts[0], 554 eight_process.use_unit_6or7ornot.disjuncts[0], 555 eight_process.use_unit_8ornot.disjuncts[0]] 556 ] 557 558 def assert_correct_disjuncts_active(nlp_model, solve_data): 559 if solve_data.master_iteration >= 1: 560 return # only checking initialization 561 iter_num = solve_data.nlp_iteration 562 disjs_should_be_active = initialize[iter_num - 1] 563 for orig_disj, soln_disj in zip( 564 solve_data.original_model.GDPopt_utils.disjunct_list, 565 nlp_model.GDPopt_utils.disjunct_list 566 ): 567 if orig_disj in disjs_should_be_active: 568 self.assertTrue(soln_disj.indicator_var.value == 1) 569 570 SolverFactory('gdpopt').solve( 571 eight_process, strategy='RIC', init_strategy='custom_disjuncts', 572 custom_init_disjuncts=initialize, 573 mip_solver=mip_solver, 574 nlp_solver=nlp_solver, 575 call_after_subproblem_feasible=assert_correct_disjuncts_active) 576 577 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 578 579 580@unittest.skipIf(not GLOA_solvers_available, 581 "Required subsolvers %s are not available" 582 % (GLOA_solvers,)) 583@unittest.skipIf(not mcpp_available(), "MC++ is not available") 584class TestGLOA(unittest.TestCase): 585 """Tests for global logic-based outer approximation.""" 586 587 def test_GDP_integer_vars(self): 588 m = ConcreteModel() 589 m.x = Var(bounds=(0, 10)) 590 m.y = Var(domain=Integers, bounds=(0, 5)) 591 m.d = Disjunction(expr=[[ 592 m.x >= m.y, m.y >= 3.5 593 ], 594 [ 595 m.x >= m.y, m.y >= 2.5 596 ]]) 597 m.o = Objective(expr=m.x) 598 SolverFactory('gdpopt').solve( 599 m, strategy='GLOA', 600 mip_solver=mip_solver, 601 nlp_solver=nlp_solver, 602 minlp_solver=minlp_solver 603 ) 604 self.assertAlmostEqual(value(m.o.expr), 3) 605 606 def test_GDP_integer_vars_infeasible(self): 607 m = ConcreteModel() 608 m.x = Var(bounds=(0, 1)) 609 m.y = Var(domain=Integers, bounds=(0, 5)) 610 m.d = Disjunction(expr=[[ 611 m.x >= m.y, m.y >= 3.5 612 ], 613 [ 614 m.x >= m.y, m.y >= 2.5 615 ]]) 616 m.o = Objective(expr=m.x) 617 res = SolverFactory('gdpopt').solve( 618 m, strategy='GLOA', 619 mip_solver=mip_solver, 620 nlp_solver=nlp_solver, 621 minlp_solver=minlp_solver 622 ) 623 self.assertEqual(res.solver.termination_condition, TerminationCondition.infeasible) 624 625 @unittest.skipUnless(license_available, "Global NLP solver license not available.") 626 def test_GLOA_8PP(self): 627 """Test the global logic-based outer approximation algorithm.""" 628 exfile = import_file( 629 join(exdir, 'eight_process', 'eight_proc_model.py')) 630 eight_process = exfile.build_eight_process_flowsheet() 631 SolverFactory('gdpopt').solve( 632 eight_process, strategy='GLOA', tee=False, 633 mip_solver=mip_solver, 634 nlp_solver=global_nlp_solver, 635 nlp_solver_args=global_nlp_solver_args 636 ) 637 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 638 639 @unittest.skipUnless(license_available, "Global NLP solver license not available.") 640 def test_GLOA_8PP_force_NLP(self): 641 """Test the global logic-based outer approximation algorithm.""" 642 exfile = import_file( 643 join(exdir, 'eight_process', 'eight_proc_model.py')) 644 eight_process = exfile.build_eight_process_flowsheet() 645 SolverFactory('gdpopt').solve( 646 eight_process, strategy='GLOA', tee=False, 647 mip_solver=mip_solver, 648 nlp_solver=global_nlp_solver, 649 nlp_solver_args=global_nlp_solver_args, 650 force_subproblem_nlp=True 651 ) 652 self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1E-2) 653 654 @unittest.skipUnless(license_available, "Global NLP solver license not available.") 655 def test_GLOA_strip_pack_default_init(self): 656 """Test logic-based outer approximation with strip packing.""" 657 exfile = import_file( 658 join(exdir, 'strip_packing', 'strip_packing_concrete.py')) 659 strip_pack = exfile.build_rect_strip_packing_model() 660 SolverFactory('gdpopt').solve( 661 strip_pack, strategy='GLOA', 662 mip_solver=mip_solver, 663 nlp_solver=global_nlp_solver, 664 nlp_solver_args=global_nlp_solver_args) 665 self.assertTrue( 666 fabs(value(strip_pack.total_length.expr) - 11) <= 1E-2) 667 668 @unittest.skipUnless(license_available, "Global NLP solver license not available.") 669 @unittest.category('expensive') 670 def test_GLOA_constrained_layout_default_init(self): 671 """Test LOA with constrained layout.""" 672 exfile = import_file( 673 join(exdir, 'constrained_layout', 'cons_layout_model.py')) 674 cons_layout = exfile.build_constrained_layout_model() 675 SolverFactory('gdpopt').solve( 676 cons_layout, strategy='GLOA', 677 mip_solver=mip_solver, 678 iterlim=36, 679 nlp_solver=global_nlp_solver, 680 nlp_solver_args=global_nlp_solver_args, 681 tee=False) 682 objective_value = value(cons_layout.min_dist_cost.expr) 683 self.assertTrue( 684 fabs(objective_value - 41573) <= 200, 685 "Objective value of %s instead of 41573" % objective_value) 686 687 def test_GLOA_ex_633_trespalacios(self): 688 """Test LOA with Francisco thesis example.""" 689 exfile = import_file(join(exdir, 'small_lit', 'ex_633_trespalacios.py')) 690 model = exfile.build_simple_nonconvex_gdp() 691 SolverFactory('gdpopt').solve( 692 model, strategy='GLOA', 693 mip_solver=mip_solver, 694 nlp_solver=global_nlp_solver, 695 nlp_solver_args=global_nlp_solver_args, 696 tee=False) 697 objective_value = value(model.obj.expr) 698 self.assertAlmostEqual(objective_value, 4.46, 2) 699 700 def test_GLOA_nonconvex_HENS(self): 701 exfile = import_file(join(exdir, 'small_lit', 'nonconvex_HEN.py')) 702 model = exfile.build_gdp_model() 703 SolverFactory('gdpopt').solve( 704 model, strategy='GLOA', 705 mip_solver=mip_solver, 706 nlp_solver=global_nlp_solver, 707 nlp_solver_args=global_nlp_solver_args, 708 tee=False) 709 objective_value = value(model.objective.expr) 710 self.assertAlmostEqual(objective_value * 1E-5, 1.14385, 2) 711 712 def test_GLOA_disjunctive_bounds(self): 713 exfile = import_file(join(exdir, 'small_lit', 'nonconvex_HEN.py')) 714 model = exfile.build_gdp_model() 715 SolverFactory('gdpopt').solve( 716 model, strategy='GLOA', 717 mip_solver=mip_solver, 718 nlp_solver=global_nlp_solver, 719 nlp_solver_args=global_nlp_solver_args, 720 calc_disjunctive_bounds=True, 721 tee=False) 722 objective_value = value(model.objective.expr) 723 self.assertAlmostEqual(objective_value * 1E-5, 1.14385, 2) 724 725 726if __name__ == '__main__': 727 unittest.main() 728