|
| 1 | +# Copyright 2017 Novo Nordisk Foundation Center for Biosustainability, |
| 2 | +# Technical University of Denmark. |
| 3 | +# |
| 4 | +# Licensed under the Apache License, Version 2.0 (the "License"); |
| 5 | +# you may not use this file except in compliance with the License. |
| 6 | +# You may obtain a copy of the License at |
| 7 | +# |
| 8 | +# http://www.apache.org/licenses/LICENSE-2.0 |
| 9 | +# |
| 10 | +# Unless required by applicable law or agreed to in writing, software |
| 11 | +# distributed under the License is distributed on an "AS IS" BASIS, |
| 12 | +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. |
| 13 | +# See the License for the specific language governing permissions and |
| 14 | +# limitations under the License. |
| 15 | + |
| 16 | + |
| 17 | +""" |
| 18 | +Interface for the GNU Linear Programming Kit (GLPK) |
| 19 | +
|
| 20 | +GLPK is an open source LP solver, with MILP capabilities. This interface exposes its GLPK's exact solver. |
| 21 | +To use GLPK you need to install the 'swiglpk' python package (with pip or from http://github.com/biosustain/swiglpk) |
| 22 | +and make sure that 'import swiglpk' runs without error. |
| 23 | +""" |
| 24 | + |
| 25 | +import logging |
| 26 | + |
| 27 | +import six |
| 28 | + |
| 29 | +from optlang.util import inheritdocstring |
| 30 | +from optlang import interface |
| 31 | +from optlang import glpk_interface |
| 32 | +from optlang.glpk_interface import _GLPK_STATUS_TO_STATUS |
| 33 | + |
| 34 | +log = logging.getLogger(__name__) |
| 35 | + |
| 36 | +from swiglpk import glp_exact, glp_create_prob, glp_get_status, \ |
| 37 | + GLP_SF_AUTO, GLP_ETMLIM, glp_adv_basis, glp_read_lp, glp_scale_prob |
| 38 | + |
| 39 | + |
| 40 | +@six.add_metaclass(inheritdocstring) |
| 41 | +class Variable(glpk_interface.Variable): |
| 42 | + def __init__(self, name, index=None, type="continuous", **kwargs): |
| 43 | + if type in ("integer", "binary"): |
| 44 | + raise ValueError("The GLPK exact solver does not support integer and mixed integer problems") |
| 45 | + super(Variable, self).__init__(name, index, type=type, **kwargs) |
| 46 | + |
| 47 | + @glpk_interface.Variable.type.setter |
| 48 | + def type(self, value): |
| 49 | + if value in ("integer", "binary"): |
| 50 | + raise ValueError("The GLPK exact solver does not support integer and mixed integer problems") |
| 51 | + super(Variable, Variable).type.fset(self, value) |
| 52 | + |
| 53 | + |
| 54 | +@six.add_metaclass(inheritdocstring) |
| 55 | +class Constraint(glpk_interface.Constraint): |
| 56 | + pass |
| 57 | + |
| 58 | + |
| 59 | +@six.add_metaclass(inheritdocstring) |
| 60 | +class Objective(glpk_interface.Objective): |
| 61 | + pass |
| 62 | + |
| 63 | + |
| 64 | +@six.add_metaclass(inheritdocstring) |
| 65 | +class Configuration(glpk_interface.Configuration): |
| 66 | + pass |
| 67 | + |
| 68 | + |
| 69 | +@six.add_metaclass(inheritdocstring) |
| 70 | +class Model(glpk_interface.Model): |
| 71 | + def _run_glp_exact(self): |
| 72 | + return_value = glp_exact(self.problem, self.configuration._smcp) |
| 73 | + glpk_status = glp_get_status(self.problem) |
| 74 | + if return_value == 0: |
| 75 | + status = _GLPK_STATUS_TO_STATUS[glpk_status] |
| 76 | + elif return_value == GLP_ETMLIM: |
| 77 | + status = interface.TIME_LIMIT |
| 78 | + else: |
| 79 | + status = _GLPK_STATUS_TO_STATUS[glpk_status] |
| 80 | + if status == interface.UNDEFINED: |
| 81 | + log.debug("Status undefined. GLPK status code returned by glp_simplex was %d" % return_value) |
| 82 | + return status |
| 83 | + |
| 84 | + def _optimize(self): |
| 85 | + # Solving inexact first per GLPK manual |
| 86 | + # Computations in exact arithmetic are very time consuming, so solving LP |
| 87 | + # problem with the routine glp_exact from the very beginning is not a good |
| 88 | + # idea. It is much better at first to find an optimal basis with the routine |
| 89 | + # glp_simplex and only then to call glp_exact, in which case only a few |
| 90 | + # simplex iterations need to be performed in exact arithmetic. |
| 91 | + status = super(Model, self)._optimize() |
| 92 | + if status != interface.OPTIMAL: |
| 93 | + return status |
| 94 | + else: |
| 95 | + status = self._run_glp_exact() |
| 96 | + |
| 97 | + if status == interface.UNDEFINED and self.configuration.presolve is True: |
| 98 | + # If presolve is on, status will be undefined if not optimal |
| 99 | + self.configuration.presolve = False |
| 100 | + status = self._run_glp_exact() |
| 101 | + self.configuration.presolve = True |
| 102 | + return status |
| 103 | + |
| 104 | + |
| 105 | +if __name__ == '__main__': |
| 106 | + import pickle |
| 107 | + |
| 108 | + x1 = Variable('x1', lb=0) |
| 109 | + x2 = Variable('x2', lb=0) |
| 110 | + x3 = Variable('x3', lb=0, ub=1, type='binary') |
| 111 | + c1 = Constraint(x1 + x2 + x3, lb=-100, ub=100, name='c1') |
| 112 | + c2 = Constraint(10 * x1 + 4 * x2 + 5 * x3, ub=600, name='c2') |
| 113 | + c3 = Constraint(2 * x1 + 2 * x2 + 6 * x3, ub=300, name='c3') |
| 114 | + obj = Objective(10 * x1 + 6 * x2 + 4 * x3, direction='max') |
| 115 | + model = Model(name='Simple model') |
| 116 | + model.objective = obj |
| 117 | + model.add([c1, c2, c3]) |
| 118 | + model.configuration.verbosity = 3 |
| 119 | + status = model.optimize() |
| 120 | + print("status:", model.status) |
| 121 | + print("objective value:", model.objective.value) |
| 122 | + |
| 123 | + for var_name, var in model.variables.items(): |
| 124 | + print(var_name, "=", var.primal) |
| 125 | + |
| 126 | + print(model) |
| 127 | + |
| 128 | + problem = glp_create_prob() |
| 129 | + glp_read_lp(problem, None, "tests/data/model.lp") |
| 130 | + |
| 131 | + solver = Model(problem=problem) |
| 132 | + print(solver.optimize()) |
| 133 | + print(solver.objective) |
| 134 | + |
| 135 | + import time |
| 136 | + |
| 137 | + t1 = time.time() |
| 138 | + print("pickling") |
| 139 | + pickle_string = pickle.dumps(solver) |
| 140 | + resurrected_solver = pickle.loads(pickle_string) |
| 141 | + t2 = time.time() |
| 142 | + print("Execution time: %s" % (t2 - t1)) |
| 143 | + |
| 144 | + resurrected_solver.optimize() |
| 145 | + print("Halelujah!", resurrected_solver.objective.value) |
0 commit comments