1# Copyright 2011 Hakan Kjellerstrand hakank@gmail.com 2# 3# Licensed under the Apache License, Version 2.0 (the "License"); 4# you may not use this file except in compliance with the License. 5# You may obtain a copy of the License at 6# 7# http://www.apache.org/licenses/LICENSE-2.0 8# 9# Unless required by applicable law or agreed to in writing, software 10# distributed under the License is distributed on an "AS IS" BASIS, 11# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 12# See the License for the specific language governing permissions and 13# limitations under the License. 14""" 15 16 3 jugs problem using MIP in Google or-tools. 17 18 A.k.a. water jugs problem. 19 20 Problem from Taha 'Introduction to Operations Research', 21 page 245f . 22 23 Compare with the CP model: 24 http://www.hakank.org/google_or_tools/3_jugs_regular 25 26 This model was created by Hakan Kjellerstrand (hakank@gmail.com) 27 Also see my other Google CP Solver models: 28 http://www.hakank.org/google_or_tools/ 29""" 30import sys 31from ortools.linear_solver import pywraplp 32 33 34def main(sol='CBC'): 35 36 # Create the solver. 37 38 print('Solver: ', sol) 39 40 # using GLPK 41 if sol == 'GLPK': 42 solver = pywraplp.Solver('CoinsGridGLPK', 43 pywraplp.Solver.GLPK_MIXED_INTEGER_PROGRAMMING) 44 else: 45 # Using CBC 46 solver = pywraplp.Solver('CoinsGridCBC', 47 pywraplp.Solver.CBC_MIXED_INTEGER_PROGRAMMING) 48 49 # 50 # data 51 # 52 n = 15 53 start = 0 # start node 54 end = 14 # end node 55 M = 999 # a large number 56 57 nodes = [ 58 '8,0,0', # start 59 '5,0,3', 60 '5,3,0', 61 '2,3,3', 62 '2,5,1', 63 '7,0,1', 64 '7,1,0', 65 '4,1,3', 66 '3,5,0', 67 '3,2,3', 68 '6,2,0', 69 '6,0,2', 70 '1,5,2', 71 '1,4,3', 72 '4,4,0' # goal! 73 ] 74 75 # distance 76 d = [[M, 1, M, M, M, M, M, M, 1, M, M, M, M, M, M], 77 [M, M, 1, M, M, M, M, M, M, M, M, M, M, M, M], 78 [M, M, M, 1, M, M, M, M, 1, M, M, M, M, M, M], 79 [M, M, M, M, 1, M, M, M, M, M, M, M, M, M, M], 80 [M, M, M, M, M, 1, M, M, 1, M, M, M, M, M, M], 81 [M, M, M, M, M, M, 1, M, M, M, M, M, M, M, M], 82 [M, M, M, M, M, M, M, 1, 1, M, M, M, M, M, M], 83 [M, M, M, M, M, M, M, M, M, M, M, M, M, M, 1], 84 [M, M, M, M, M, M, M, M, M, 1, M, M, M, M, M], 85 [M, 1, M, M, M, M, M, M, M, M, 1, M, M, M, M], 86 [M, M, M, M, M, M, M, M, M, M, M, 1, M, M, M], 87 [M, 1, M, M, M, M, M, M, M, M, M, M, 1, M, M], 88 [M, M, M, M, M, M, M, M, M, M, M, M, M, 1, M], 89 [M, 1, M, M, M, M, M, M, M, M, M, M, M, M, 1], 90 [M, M, M, M, M, M, M, M, M, M, M, M, M, M, M]] 91 92 # 93 # variables 94 # 95 96 # requirements (right hand statement) 97 rhs = [solver.IntVar(-1, 1, 'rhs[%i]' % i) for i in range(n)] 98 99 x = {} 100 for i in range(n): 101 for j in range(n): 102 x[i, j] = solver.IntVar(0, 1, 'x[%i,%i]' % (i, j)) 103 104 out_flow = [solver.IntVar(0, 1, 'out_flow[%i]' % i) for i in range(n)] 105 in_flow = [solver.IntVar(0, 1, 'in_flow[%i]' % i) for i in range(n)] 106 107 # length of path, to be minimized 108 z = solver.Sum( 109 [d[i][j] * x[i, j] for i in range(n) for j in range(n) if d[i][j] < M]) 110 111 # 112 # constraints 113 # 114 115 for i in range(n): 116 if i == start: 117 solver.Add(rhs[i] == 1) 118 elif i == end: 119 solver.Add(rhs[i] == -1) 120 else: 121 solver.Add(rhs[i] == 0) 122 123 # outflow constraint 124 for i in range(n): 125 solver.Add( 126 out_flow[i] == solver.Sum([x[i, j] for j in range(n) if d[i][j] < M])) 127 128 # inflow constraint 129 for j in range(n): 130 solver.Add( 131 in_flow[j] == solver.Sum([x[i, j] for i in range(n) if d[i][j] < M])) 132 133 # inflow = outflow 134 for i in range(n): 135 solver.Add(out_flow[i] - in_flow[i] == rhs[i]) 136 137 # objective 138 objective = solver.Minimize(z) 139 140 # 141 # solution and search 142 # 143 solver.Solve() 144 145 print() 146 print('z: ', int(solver.Objective().Value())) 147 148 t = start 149 while t != end: 150 print(nodes[t], '->', end=' ') 151 for j in range(n): 152 if x[t, j].SolutionValue() == 1: 153 print(nodes[j]) 154 t = j 155 break 156 157 print() 158 print('walltime :', solver.WallTime(), 'ms') 159 if sol == 'CBC': 160 print('iterations:', solver.Iterations()) 161 162 163if __name__ == '__main__': 164 165 sol = 'CBC' 166 if len(sys.argv) > 1: 167 sol = sys.argv[1] 168 if sol != 'GLPK' and sol != 'CBC': 169 print('Solver must be either GLPK or CBC') 170 sys.exit(1) 171 172 main(sol) 173