1111
1212import sys
1313
14- # Workaround for issues with python not being installed as a framework on mac
15- # by using a different backend.
16- if sys .platform == "darwin" : # OS X
17- import matplotlib as mpl
18- mpl .use ('Agg' )
19- del mpl
20-
21- import matplotlib .pyplot as plt
2214from math import sqrt , log
2315from itertools import product
2416from mip import Model , xsum , minimize , OptimizationStatus
3628C = [1 , 2 , 3 , 4 , 5 , 6 , 7 , 8 , 9 , 10 ]
3729
3830# position of clients
39- pc = {1 : (94 , 10 ), 2 : (57 , 26 ), 3 : (74 , 44 ), 4 : (27 , 51 ), 5 : (78 , 30 ), 6 : (23 , 30 ),
40- 7 : (20 , 72 ), 8 : (3 , 27 ), 9 : (5 , 39 ), 10 : (51 , 1 )}
31+ pc = {
32+ 1 : (94 , 10 ),
33+ 2 : (57 , 26 ),
34+ 3 : (74 , 44 ),
35+ 4 : (27 , 51 ),
36+ 5 : (78 , 30 ),
37+ 6 : (23 , 30 ),
38+ 7 : (20 , 72 ),
39+ 8 : (3 , 27 ),
40+ 9 : (5 , 39 ),
41+ 10 : (51 , 1 ),
42+ }
4143
4244# demands
4345d = {1 : 302 , 2 : 273 , 3 : 275 , 4 : 266 , 5 : 287 , 6 : 296 , 7 : 297 , 8 : 310 , 9 : 302 , 10 : 309 }
4446
45- # plotting possible plant locations
46- for i , p in pf .items ():
47- plt .scatter ((p [0 ]), (p [1 ]), marker = "^" , color = "purple" , s = 50 )
48- plt .text ((p [0 ]), (p [1 ]), "$f_%d$" % i )
49-
50- # plotting location of clients
51- for i , p in pc .items ():
52- plt .scatter ((p [0 ]), (p [1 ]), marker = "o" , color = "black" , s = 15 )
53- plt .text ((p [0 ]), (p [1 ]), "$c_{%d}$" % i )
54-
55- plt .text ((20 ), (78 ), "Region 1" )
56- plt .text ((70 ), (78 ), "Region 2" )
57- plt .plot ((50 , 50 ), (0 , 80 ))
58-
59- dist = {(f , c ): round (sqrt ((pf [f ][0 ] - pc [c ][0 ]) ** 2 + (pf [f ][1 ] - pc [c ][1 ]) ** 2 ), 1 )
60- for (f , c ) in product (F , C ) }
47+ dist = {
48+ (f , c ): round (sqrt ((pf [f ][0 ] - pc [c ][0 ]) ** 2 + (pf [f ][1 ] - pc [c ][1 ]) ** 2 ), 1 )
49+ for (f , c ) in product (F , C )
50+ }
6151
6252m = Model ()
6353
8272 D = 6 # nr. of discretization points, increase for more precision
8373 v = [c [f ] * (v / (D - 1 )) for v in range (D )] # points
8474 # non-linear function values for points in v
85- vn = [0 if k == 0 else 1520 * log (v [k ]) for k in range (D )]
75+ vn = [0 if k == 0 else 1520 * log (v [k ]) for k in range (D )]
8676 # w variables
8777 w = [m .add_var () for v in range (D )]
8878 m += xsum (w ) == 1 # convexification
9888
9989# objective function
10090m .objective = minimize (
101- xsum (dist [i , j ] * x [i , j ] for (i , j ) in product (F , C )) + xsum (y [i ] for i in F ) )
91+ xsum (dist [i , j ] * x [i , j ] for (i , j ) in product (F , C )) + xsum (y [i ] for i in F )
92+ )
10293
10394m .optimize ()
10495
105- plt .savefig ("location.pdf" )
106-
10796if m .num_solutions :
10897 print ("Solution with cost {} found." .format (m .objective_value ))
10998 print ("Facilities capacities: {} " .format ([z [f ].x for f in F ]))
11099 print ("Facilities cost: {}" .format ([y [f ].x for f in F ]))
111100
112- # plotting allocations
113- for (i , j ) in [(i , j ) for (i , j ) in product (F , C ) if x [(i , j )].x >= 1e-6 ]:
114- plt .plot (
115- (pf [i ][0 ], pc [j ][0 ]), (pf [i ][1 ], pc [j ][1 ]), linestyle = "--" , color = "darkgray"
116- )
117-
118- plt .savefig ("location-sol.pdf" )
119-
120101# sanity checks
121102opt = 99733.94905406
122103if m .status == OptimizationStatus .OPTIMAL :
123104 assert abs (m .objective_value - opt ) <= 0.01
124105elif m .status == OptimizationStatus .FEASIBLE :
125106 assert m .objective_value >= opt - 0.01
126107else :
127- assert m .status not in [OptimizationStatus .INFEASIBLE , OptimizationStatus .UNBOUNDED ]
108+ assert m .status not in [OptimizationStatus .INFEASIBLE , OptimizationStatus .UNBOUNDED ]
0 commit comments