Added CPLEX solver model in MIP class

Week 6 - Facility Location/facility/EnumSettings.py

@@ -25,4 +25,8 @@ class SolvingParadigm(Enum):
 class InitialSolutionFunction(Enum):
     Radius = "Radius"
     Euclidean = "Euclidean"
-    Manhatan = "Manhatan"
+    Manhatan = "Manhatan"
+
+class MipSolver(Enum):
+    SCIP = "SCIP"
+    CPLEX = "CPLEX"

Week 6 - Facility Location/facility/LNS.py

@@ -1,7 +1,7 @@
 from Forest import Forest
 from Util import Util
 from MIP import MIP
-from EnumSettings import Strategy,ImprovementType,SolvingParadigm,InitialSolutionFunction
+from EnumSettings import Strategy,ImprovementType,SolvingParadigm,InitialSolutionFunction,MipSolver
 from Preprocessing import Preprocessing
 from Tree import Tree
 import math
@@ -24,7 +24,7 @@ def __init__(self,facilities,customers,params):
         self.customers = customers
         self.subproblemSolutionForest = Forest()
         self.currentIteration = 0
-        self.mip = MIP(facilities,customers,"Instance_%s_%s" %(len(facilities),len(customers)))
+        self.mip = MIP(facilities,customers,"Instance_%s_%s" %(len(facilities),len(customers)),params["mipSolver"])
         self.facilitiesCount = len(facilities)
         self.quantiles = []
         self.params = params
@@ -109,8 +109,8 @@ def __destroy(self,cluster):
                 for node in self.subproblemSolutionForest.getTrees().get(facilityIndex).getNodes().values():
                     clusterDemand = clusterDemand + node.demand
 
-        candidateFacilities = self.__getCandidateFacilities(cluster,clusterDemand,Util.truncate(self.quantiles[self.currentIteration],5))
-
+        #candidateFacilities = self.__getCandidateFacilities(cluster,clusterDemand,Util.truncate(self.quantiles[self.currentIteration],5))
+        candidateFacilities = copy.deepcopy(cluster)
         print("Current Forest: %s/%s - Candidate Facilities: %s/%s"%(self.subproblemSolutionForest.getTreesCount(),self.subproblemSolutionForest.getTotalNodes(),len(candidateFacilities),len(cluster)))
 
         reassignmentCandidates = Forest()
@@ -138,8 +138,7 @@ def __repair(self,candidatesFacility,candidatesCustomer):
             print("=============================")
             print("Repair Method Started...")
         self.mip.clear()
-        self.mip.initialize(candidatesFacility,candidatesCustomer,"Instance_%s_%s" %(len(candidatesFacility),len(candidatesCustomer)))
-        self.mip.createModel()
+        self.mip.initialize(candidatesFacility,candidatesCustomer,"Instance_%s_%s" %(len(candidatesFacility),len(candidatesCustomer)),self.params["mipSolver"])
         obj,assignments,status = self.mip.optimize(self.params["mipTimeLimit"])
 
         if(self.DEBUG_MESSAGES):

Week 6 - Facility Location/facility/MIP.py

@@ -1,12 +1,16 @@
 from pyscipopt import Model, quicksum
 import math
 from Preprocessing import Preprocessing
+from EnumSettings import MipSolver
+import docplex.mp.model as cpx
+from docplex.util.status import JobSolveStatus
+
 
 class MIP:
-    DEBUG_MESSAGES = False
+    DEBUG_MESSAGES = True
 
-    def __init__(self, f, c, instanceName): 
-        self.initialize(f,c,instanceName)
+    def __init__(self, f, c, instanceName,mipSolver): 
+        self.initialize(f,c,instanceName,mipSolver)
 
     def clear(self):
         self.facilities = [] 
@@ -15,15 +19,25 @@ def clear(self):
         self.varFacilityAssignment = {}
         self.varCustomerAssignment = {}
 
-    def initialize(self, f, c, instanceName):
+    def initialize(self, f, c, instanceName,mipSolver):
         self.clear()
         self.facilities = f 
         self.customers = c 
         self.instanceName = instanceName
         self.varFacilityAssignment = {}
         self.varCustomerAssignment = {}
+        self.solver = mipSolver
+
+    def optimize(self,timeLimit):       
+        if(self.solver == MipSolver.SCIP):
+            self.__createModelSCIP()
+            return self.__optimizeSCIP(timeLimit)
+
+        elif (self.solver == MipSolver.CPLEX):
+            self.__createModelCPLEX()
+            return self.__optimizeCPLEX(timeLimit)
 
-    def createModel(self):
+    def __createModelSCIP(self):
         self.model = Model(self.instanceName)
         print("Creating Variables...")
         #Variables
@@ -54,9 +68,71 @@ def createModel(self):
         self.model.setObjective(quicksum(self.varFacilityAssignment[facility.index]*facility.setup_cost for facility in self.facilities) + quicksum(Preprocessing.getEuclideanDistance(facility.location,customer.location)*self.varCustomerAssignment[facility.index,customer.index] for facility in self.facilities for customer in self.customers),"minimize")
         self.model.data = self.varFacilityAssignment, self.varCustomerAssignment
 
-    def optimize(self,timeLimit):
+    def __createModelCPLEX(self):
+        self.model = cpx.Model(self.instanceName)
+        print("Creating Variables...")
+        #Variables
+        for f in self.facilities:
+            self.varFacilityAssignment[f.index] = self.model.binary_var(name="facility-%s" % f.index)
+            for c in self.customers:
+                #Demand is binary because each customer must be served by exaclty one facility
+                self.varCustomerAssignment[f.index,c.index] = self.model.binary_var(name="demand-(%s,%s)" % (f.index,c.index))
+
+        print("Creating Constraints...")
+        #Constraints
+        #Ensure all customers are assigned to one facility
+        for customer in self.customers:
+            self.model.add_constraint(ct=self.model.sum(self.varCustomerAssignment[facility.index,customer.index] for facility in self.facilities) == 1,ctname="Demand(%s)"% customer.index)
+
+        #Ensure the demand carried by the facility is at most its capacity
+        for facility in self.facilities:
+              self.model.add_constraint(ct=self.model.sum(self.varCustomerAssignment[facility.index,customer.index]*customer.demand for customer in self.customers) <= facility.capacity*self.varFacilityAssignment[facility.index],ctname="Capacity(%s)" % facility.index)
+
+        #Strong Formulation
+        for facility in self.facilities:
+            for customer in self.customers:
+                self.model.add_constraint(ct=self.varCustomerAssignment[facility.index,customer.index] <=  facility.capacity*self.varFacilityAssignment[facility.index],ctname="Strong(%s,%s)"%(facility.index,customer.index))
+                self.model.add_constraint(self.varCustomerAssignment[facility.index,customer.index] >= 0,ctname="Strong2(%s,%s)"%(facility.index,customer.index))
+
+        print("Creating Objective Function...")
+        #Objective Function
+        objective = self.model.sum(self.varFacilityAssignment[facility.index]*facility.setup_cost for facility in self.facilities) + self.model.sum(Preprocessing.getEuclideanDistance(facility.location,customer.location)*self.varCustomerAssignment[facility.index,customer.index] for facility in self.facilities for customer in self.customers)
+        self.model.minimize(objective)
+
+
+    def __optimizeCPLEX(self,timeLimit):
         print("Instace: %s" % self.instanceName)
+        if not self.DEBUG_MESSAGES:
+            self.model.print_information()
+
+        self.model.parameters.timelimit = timeLimit
+        self.model.parameters.threads = 8
+
+        print("MIP - Optimizing...")
+        solution = self.model.solve(log_output=self.DEBUG_MESSAGES)
+
+        if solution is not None:
+            print("Instace: %s solved." % self.instanceName)
+        else:
+            print("Instace: %s is infeasible" % self.instanceName)
+            return None
+
 
+        assignments = []
+
+        for f in self.facilities:
+            for c in self.customers:
+                if(self.varCustomerAssignment[f.index,c.index].solution_value == 1):
+                    assignments.append((f.index,c.index))
+
+        obj = self.model.objective_value
+        status = "No optimal"
+        if self.model.get_solve_status() == JobSolveStatus.OPTIMAL_SOLUTION:
+            status = "optimal"
+        return  obj,assignments,status
+
+    def __optimizeSCIP(self,timeLimit):
+        print("Instace: %s" % self.instanceName)
         if not self.DEBUG_MESSAGES:
             self.model.hideOutput()
 

Week 6 - Facility Location/facility/ParametersConfiguration.py

@@ -1,5 +1,5 @@
 from Util import Util
-from EnumSettings import Strategy,ImprovementType,SolvingParadigm,InitialSolutionFunction
+from EnumSettings import Strategy,ImprovementType,SolvingParadigm,InitialSolutionFunction,MipSolver
 import time
 import datetime
 
@@ -53,11 +53,12 @@ def __DefaultSetup(self,instanceSize):
         self.params["improvementType"] = ImprovementType.First
         self.params["executionTimeLimit"] = Util.getTimeInSeconds(4,50,0) #4 hours and 30 minutes of time limit
         self.params["noImprovementTimeLimit"] = Util.getTimeInSeconds(0,20,0)
-        self.params["mipTimeLimit"] = Util.getTimeInSeconds(1,0,0) #30 limits for each Mip Execution
+        self.params["mipTimeLimit"] = Util.getTimeInSeconds(4,0,0) 
         self.params["strategy"] = Strategy.Default
         self.params["paradigm"] = SolvingParadigm.MIP
         self.params["quantile_intervals"] = self.__getQuantilesIntervals()
-        self.params["initialSolutionFunction"] = InitialSolutionFunction.Radius
+        self.params["initialSolutionFunction"] = InitialSolutionFunction.Euclidean
+        self.params["mipSolver"] = MipSolver.CPLEX
 
     def __AlphaSetup(self,instanceSize):
         self.__DefaultSetup(instanceSize)

Week 6 - Facility Location/facility/Preprocessing.py

@@ -84,7 +84,7 @@ def getRadiusDistanceInitialSolution(facilities,customers,clusters):
         for customer in customers:
             customersToBeAssigned[customer.index] = customer.index
 
-        quantileIntervalSize = len(facilities[0].distance_quantiles)
+        quantileIntervalSize = len(list(facilities.values())[0].distance_quantiles)
         quantileIntervalCount = 0
         factor = 1.00
         additional = 0.05

Week 6 - Facility Location/facility/solver.py

@@ -54,8 +54,7 @@ def solve_it(input_data):
     print("Instace Size: %s || Strategy: %s || Paradigm: %s || Improvement Type: %s" % (paramsConfig.instanceSize,params["strategy"],params["paradigm"],params["improvementType"]))
     print("============================================================================================================================================================")
     if(params["paradigm"] == SolvingParadigm.MIP):
-        instance = MIP(facilities,customers,"Instance_%s_%s" %(facility_count,customer_count))
-        instance.createModel()
+        instance = MIP(facilities,customers,"Instance_%s_%s" %(facility_count,customer_count),params["mipSolver"])
         obj,assignments,_ = instance.optimize(params["mipTimeLimit"])
         output_data = '%.2f' % obj + ' ' + str(1) + '\n'
         output_data += ' '.join(map(str,Util.formatSolutionFromMIP(assignments)))