First working solution with local search. Close neighborhood.

facility/facility.py

@@ -1,95 +1,172 @@
 #!/usr/bin/env python3
 
 import math
+import logging
 from collections import namedtuple
 
 Point = namedtuple("Point", ['x', 'y'])
-Facility = namedtuple("Facility", ['index', 'setup_cost', 'capacity', 'location'])
+length = None
-Customer = namedtuple("Customer", ['index', 'demand', 'location'])
+EPSILON = 0.001
+
+
+class Length(object):
+    """ Length is a helper object o get the distance between
+        customers and facilities. It also has functions
+        to return the facilities that are closest to a customer. """
+
+    def __init__(self, facilities, customers):
+        self.customer_to_facility = [[self.length(c.location, f.location)
+                                      for f in facilities]
+                                      for c in customers]
+
+        self.customer_closes_facility = [sorted([f for f in facilities],
+                                                 key=lambda f: self.get(c, f))
+                                         for c in customers]
+
+    def get(self, customer, facility):
+        c_idx = customer.index if type(customer) is Customer else customer
+        f_idx = facility.index if type(facility) is Facility else facility
+        return self.customer_to_facility[c_idx][f_idx]
+
+    def length(self, point1, point2):
+        return math.sqrt((point1.x - point2.x)**2 + (point1.y - point2.y)**2)
+
+    def get_facilities(self, customer):
+        """ Returns closest facility in increasing order of distance. """
+        return self.customer_closes_facility[customer.index]
+
+    def get_feasible_facilities(self, customer, facilities):
+        return (f for f in self.get_facilities(customer)
+                if f.remaining_capacity >= customer.demand)
+
+    def get_feasible_open_facilities(self, customer, facilities):
+        return (f for f in self.get_facilities(customer)
+                if f.remaining_capacity >= customer.demand
+                if f.is_open)
+
+
+class Facility(object):
+    def __init__(self, index, setup_cost, capacity, location):
+        self.index = index
+        self.setup_cost = setup_cost
+        self.capacity = capacity
+        self.location = location
+        self.is_open = False
+        self.remaining_capacity = capacity
+        self.customers = set()
+
+    def remove(self, customer):
+        logging.debug(f"From {self} remove {customer}.")
+        if not customer in self.customers:
+            raise ValueError(f"{customer} not connected to {self}.")
+
+        delta = 0
+        self.customers.remove(customer)
+        self.remaining_capacity += customer.demand
+        customer.facility = None
+
+        if not self.customers and self.is_open:
+            logging.debug(f"{self} is empty but open.")
+
+        delta -= length.get(customer.index, self.index)
+        return delta
+
+    def add(self, customer):
+        logging.debug(f"To {self} add {customer}.")
+        if not self.is_open:
+            raise ValueError(f"Cannot connect {customer} to not open {self}.")
+
+        if customer.demand > self.remaining_capacity:
+            raise ValueError(f"Capacity of {self} too low for {customer}")
+
+        delta = 0
+        if (other_facility := customer.facility):
+            delta += other_facility.remove(customer)
+
+        self.remaining_capacity -= customer.demand
+        self.customers.add(customer)
+        delta += length.get(customer.index, self.index)
+        customer.facility = self
+        return delta
+
+    def set_open(self):
+        logging.debug(f"Open {self}.")
+        self.is_open = True
+        return self.setup_cost
+
+    def set_not_open(self):
+        logging.debug(f"Close {self}.")
+        delta = 0
+        if self.customers:
+            raise ValueError(f"Cannot close {self} with {self.customers}.")
+        self.is_open = False
+        return -self.setup_cost
+
+    def remove_all_and_close(self):
+        logging.debug(f"{self} remove all and close.")
+        delta = 0
+        for customer in list(self.customers):
+            delta += self.remove(customer)
+        delta += self.set_not_open()
+        return delta
+
+    def __str__(self):
+        cap = f"{self.remaining_capacity}/{self.capacity}"
+        status = "O" if self.is_open else "C"
+        s = f"F({self.index}, {cap}, {status}, {self.setup_cost})"
+        return s
+
+    def __repr__(self):
+        return self.__str__()
+
+
+class Customer(object):
+    def __init__(self, index, demand, location):
+        self.index = index
+        self.demand = demand
+        self.location = location
+        self.facility = None
+
+    def __str__(self):
+        con = "C" if self.facility is not None else "NC!"
+        s = f"C({self.index}, {con})"
+        return s
+
+    def __repr__(self):
+        return self.__str__()
 
 
 class Solution(object):
 
     def __init__(self, facilities, customers):
-        self.facilities = facilities
+        self.fs = facilities
-        self.customers = customers
+        self.cs = customers
-
         self.cost = 0
-        self.facility_connected_customers = [set() for _ in facilities]
-        self.facility_remaining_capacity = [f.capacity for f in facilities]
-        self.customer_to_facility = [None for _ in customers]
 
-    def connect(self, customer_index, facility_index):
+    def validate(self):
-        customer = self.customers[customer_index]
+        cost = 0
-        facility = self.facilities[facility_index]
 
-        # If customers is already connected handle disconnect properly.
+        for c in self.cs:
-        if (connected_facitlity_index := self.customer_to_facility[customer_index]):
+            if c.facility is None:
-            self.disconnect(customer_index, connected_facitlity_index)
+                raise Exception(f"{c} not connected.")
+            cost += length.get(c.index, c.facility.index)
 
-        # If facility is currently not used we have to set it up.
+        for f in self.fs:
-        if not self.facility_connected_customers[facility_index]:
+            if f.remaining_capacity < 0:
-            self.cost += facility.setup_cost
+                raise Exception(f"{f} exceeds capacity.")
 
-        self.facility_connected_customers[facility_index].add(customer_index)
+            if f.customers and not f.is_open:
-        if self.facility_remaining_capacity[facility_index] < customer.demand:
+                raise Exception(f"{f} has customers, but is not open.")
-            raise Exception(f"Cannot connect {customer} to {facility}.")
-        self.facility_remaining_capacity[facility_index] -= customer.demand
 
-        self.customer_to_facility[customer_index] = facility_index
+            if not f.customers and f.is_open:
-        self.cost += length(facility.location, customer.location)
+                raise Exception(f"{f} has no customers, but is open.")
 
-    def disconnect(self, customer_index, facility_index):
+            if f.is_open:
-        customer = self.customers[customer_index]
+                cost += f.setup_cost
-        facility = self.facilities[facility_index]
-        self.cost -= length(facility.location, customer.location)
 
-        self.facility_connected_customers[facility_index].remove(customer_index)
+        if abs(cost - self.cost) > EPSILON:
-        self.facility_remaining_capacity[facility_index] += customer.demand
+            raise Exception(f"Running cost {self.cost} unequal to {cost}.")
 
-        self.customer_to_facility[customer_index] = None
-
-        if not self.facility_connected_customers[facility_index]:
-            self.cost -= self.facilities[facility_index].setup_cost
-
-    def get_feasible_facilities(self, customer_index):
-        customer = self.customers[customer_index]
-        facility_indices = [f.index
-                            for f in self.facilities
-                            if self.facility_remaining_capacity[f.index] >= customer.demand]
-        if not facility_indices:
-            raise Exception("No feasible facilities.")
-
-        def key(facility_index):
-            cost = 0
-            facility = self.facilities[facility_index]
-            # If there are customers yet we have to open it.
-            if not self.facility_connected_customers[facility_index]:
-                cost += facility.setup_cost
-            cost += length(customer.location, facility.location)
-            return cost
-        facility_indices.sort(key=key)
-        return facility_indices
-
-    def is_valid(self):
-        for customer in self.customers:
-            if self.customer_to_facility[customer.index] is None:
-                raise Exception(f"{customer} not connected.")
-        for facility in self.facilities:
-            if self.facility_remaining_capacity[facility.index] < 0:
-                raise Exception(f"{facility} exceeds capacity.")
-
-        cost = sum([f.setup_cost
-                    for f in self.facilities
-                    if self.facility_connected_customers[f.index]])
-
-        for customer in self.customers:
-            facility = self.facilities[self.customer_to_facility[customer.index]]
-
-            cost += length(facility.location, customer.location)
-
-        if abs(cost - self.cost) > 0.00001:
-            raise Exception(f"Running cost {self.cost} unequal to actual cost {cost}.")
         return True
 
     def plot_map(self):
@@ -101,76 +178,133 @@ class Solution(object):
 
         figure = plt.figure()
 
-        for f in self.facilities:
+        for f in self.fs:
             x, y = f.location
-            color = 'ro' if self.facility_connected_customers[f.index] else 'go'
+            color = 'ro' if f.is_open else 'go'
             plt.plot(x, y, color)
-            rem_cap = self.facility_remaining_capacity[f.index]
+            plt.text(x, y, f"{f}")
-            plt.text(x, y, f"  F({f.index}, {f.setup_cost}, {rem_cap}/{f.capacity})")
 
-        for c in self.customers:
+        for c in self.cs:
             x, y = c.location
             plt.plot(x, y, 'bx')
-            plt.text(x, y, f"  C({c.index}, {c.demand})")
+            plt.text(x, y, f"{c}")
-            if (f_index := self.customer_to_facility[c.index]) is not None:
+
-                f = self.facilities[f_index]
+            if (f := c.facility) is not None:
                 x_f, y_f = f.location
                 plt.plot([x, x_f], [y, y_f], 'b-')
         plt.show()
 
-    def get_facilities_by_customers(self):
-        facility_indices = [f.index for f in self.facilities
-                            if self.facility_connected_customers[f.index]]
-        def key(facility_index):
-            return len(self.facility_connected_customers[facility_index])
-        facility_indices.sort(key=key)
-        return facility_indices
-
     def to_output_data(self):
         # calculate the cost of the solution
-        self.is_valid()
+        self.validate()
         obj = self.cost
+        customer_to_facility = [c.facility.index for c in self.cs]
         # prepare the solution in the specified output format
         output_data = '%.2f' % obj + ' ' + str(0) + '\n'
-        output_data += ' '.join(map(str, self.customer_to_facility))
+        output_data += ' '.join(map(str, customer_to_facility))
         return output_data
 
+    def build_trivial(self):
+        facility = self.fs[0]
+        self.cost += facility.set_open()
+        for customer in self.cs:
+            if facility.remaining_capacity < customer.demand:
+                facility = self.fs[facility.index + 1]
+                self.cost += facility.set_open()
+            self.cost += facility.add(customer)
+        return self
 
-def length(point1, point2):
+    def build_greedy(self):
-    return math.sqrt((point1.x - point2.x)**2 + (point1.y - point2.y)**2)
+        def connect_to_closest_facility(customer):
+            cost = 0
+            for f in length.get_feasible_facilities(customer, self.fs):
+                if not f.is_open:
+                    cost += f.set_open()
+                cost += f.add(customer)
+                return cost
+            raise Exception("No feasible facilities for {customer}.")
+
+        for customer in self.cs:
+            self.cost += connect_to_closest_facility(customer)
+        return self
+
+    def reconnect_greedy(self, customers):
+        delta = 0
+        not_connected = []
+
+        def connect_better_facility(customer):
+            current_facility = customer.facility
+            if current_facility:
+                current_length = length.get(customer, current_facility)
+            else:
+                current_length = float("inf")
+            for f in length.get_feasible_open_facilities(customer, self.fs):
+                new_length = length.get(customer, f)
+                if new_length < current_length:
+                    logging.debug(f"{f} is better for {customer}.")
+                    return f.add(customer)
+                elif new_length > current_length:
+                    return 0
+            not_connected.append(customer)
+            return 0
+        delta += sum([connect_better_facility(c) for c in customers])
+        return delta, not_connected
 
 
-def build_trivial_solution(solution):
+    def close_facility(self, facility):
-    # build a trivial solution
+        logging.debug(f"Closing {facility}.")
-    # pack the facilities one by one until all the customers are served
+        original_cost = self.cost
-    facility_index = 0
+
-    for customer in solution.customers:
+        customers = list(facility.customers)
-        if solution.facility_remaining_capacity[facility_index] >= customer.demand:
+        self.cost += facility.remove_all_and_close()
-            solution.connect(customer.index, facility_index)
+        delta, not_connected = self.reconnect_greedy(customers)
+        self.cost += delta
+
+        if not_connected:
+            logging.info("Not all customers connected. Restore.")
+        elif self.cost < original_cost:
+            delta = original_cost - self.cost
+            logging.info(f"Close {facility} resulted in improvement {delta}.")
+            return
         else:
-            facility_index += 1
+            logging.debug(f"No improvement. Restore.")
-            solution.connect(customer.index, facility_index)
-    return solution
 
+        self.cost += facility.set_open()
+        delta, not_connected = self.reconnect_greedy(customers)
+        self.cost += delta
+        assert(not not_connected)
+        assert(abs(original_cost - self.cost) < EPSILON)
 
-def build_greedy_solution(solution):
+    def open_facility(self, facility):
-    for customer in solution.customers:
+        pass
-        facility_index = solution.get_feasible_facilities(customer.index)[0]
+
-        solution.connect(customer.index, facility_index)
+    def local_search(self):
-    return solution
+        fs = [f for f in self.fs if f.is_open]
+        fs.sort(key=lambda f: len(f.customers))
+
+        for f in fs:
+            self.close_facility(f)
+
+        #fs = [f for f in self.fs if not f.is_open]
+        #fs.sort(key=lambda f: f.setup_cost)
+        #for
 
 
 def solve_it(input_data):
+    global length
     facilities, customers = parse(input_data)
+    length = Length(facilities, customers)
     solution = Solution(facilities, customers)
-    build_greedy_solution(solution)
+    solution.build_greedy()
-    solution.plot_map()
+    solution.reconnect_greedy(solution.cs)
+    solution.local_search()
+    # solution.plot_map()
     output_data = solution.to_output_data()
     return output_data
 
 
 def main():
-    file_location = "data/fl_3_1"
+    file_location = "data/fl_100_7"
     with open(file_location, 'r') as input_data_file:
         input_data = input_data_file.read()
     print(solve_it(input_data))
@@ -200,5 +334,6 @@ def parse(input_data):
 
 
 if __name__ == "__main__":
+    logging.basicConfig(level=logging.INFO, format='%(message)s')
     main()