diff --git a/benchmarks/run.py b/benchmarks/run.py index 811fdd7..907c2b6 100644 --- a/benchmarks/run.py +++ b/benchmarks/run.py @@ -4,7 +4,9 @@ from multiprocessing import Pool, cpu_count from pathlib import Path from typing import List, Dict, Union + from networkx import DiGraph + from benchmarks.augerat_dataset import AugeratDataSet from benchmarks.solomon_dataset import SolomonDataSet from benchmarks.utils.csv_table import CsvTable @@ -202,7 +204,7 @@ def _run_single_problem(path_to_instance: Path, **kwargs): def main(): - """ Run parallel or series""" + """Run parallel or series""" if SERIES: run_series() else: diff --git a/benchmarks/tests/graph_issue101 b/benchmarks/tests/graph_issue101 new file mode 100644 index 0000000..de81fb8 Binary files /dev/null and b/benchmarks/tests/graph_issue101 differ diff --git a/benchmarks/tests/pytest.ini b/benchmarks/tests/pytest.ini index 5e0b9f0..ad1f990 100644 --- a/benchmarks/tests/pytest.ini +++ b/benchmarks/tests/pytest.ini @@ -1,5 +1,9 @@ [pytest] +log_cli = 1 +log_cli_level = INFO +log_cli_format = %(asctime)s [%(levelname)8s] %(message)s (%(filename)s:%(lineno)s) +log_cli_date_format=%Y-%m-%d %H:%M:%S python_files = test_*.py #testpaths = tests/ filterwarnings = - ignore::DeprecationWarning \ No newline at end of file + ignore::DeprecationWarning diff --git a/benchmarks/tests/test_cvrptw_solomon.py b/benchmarks/tests/test_cvrptw_solomon.py index 25eee97..5c05acd 100644 --- a/benchmarks/tests/test_cvrptw_solomon.py +++ b/benchmarks/tests/test_cvrptw_solomon.py @@ -4,19 +4,18 @@ class TestsSolomon: - def setup(self): """ Solomon instance c101, 25 first nodes only including depot """ - data = SolomonDataSet(path="benchmarks/data/cvrptw/", - instance_name="C101.txt", - n_vertices=25) + data = SolomonDataSet( + path="benchmarks/data/cvrptw/", instance_name="C101.txt", n_vertices=25 + ) self.G = data.G self.n_vertices = 25 - self.prob = VehicleRoutingProblem(self.G, - load_capacity=data.max_load, - time_windows=True) + self.prob = VehicleRoutingProblem( + self.G, load_capacity=data.max_load, time_windows=True + ) initial_routes = [ ["Source", 13, 17, 18, 19, 15, 16, 14, 12, 1, "Sink"], ["Source", 20, 24, 23, 22, 21, "Sink"], @@ -26,7 +25,7 @@ def setup(self): # Set repeating solver arguments self.solver_args = { "pricing_strategy": "BestPaths", - "initial_routes": initial_routes + "initial_routes": initial_routes, } def test_setup_instance_name(self): @@ -40,43 +39,18 @@ def test_setup_nodes(self): assert len(self.G.nodes()) == self.n_vertices + 1 def test_setup_edges(self): - assert len( - self.G.edges()) == self.n_vertices * (self.n_vertices - 1) + 1 + assert len(self.G.edges()) == self.n_vertices * (self.n_vertices - 1) + 1 def test_subproblem_lp(self): # benchmark result # e.g., in Feillet et al. (2004) self.prob.solve(**self.solver_args, cspy=False) assert round(self.prob.best_value, -1) in [190, 200] - - def test_schedule_lp(self): - 'Tests whether final schedule is time-window feasible' - self.prob.solve(**self.solver_args, cspy=False) - # Check arrival times - for k1, v1 in self.prob.arrival_time.items(): - for k2, v2 in v1.items(): - assert (self.G.nodes[k2]["lower"] <= v2) - assert (v2 <= self.G.nodes[k2]["upper"]) - # Check departure times - for k1, v1 in self.prob.departure_time.items(): - for k2, v2 in v1.items(): - assert (self.G.nodes[k2]["lower"] <= v2) - assert (v2 <= self.G.nodes[k2]["upper"]) + self.prob.check_arrival_time() + self.prob.check_departure_time() def test_subproblem_cspy(self): self.prob.solve(**self.solver_args, cspy=True) assert round(self.prob.best_value, -1) in [190, 200] - - def test_schedule_cspy(self): - 'Tests whether final schedule is time-window feasible' - self.prob.solve(**self.solver_args) - # Check departure times - for k1, v1 in self.prob.departure_time.items(): - for k2, v2 in v1.items(): - assert (self.G.nodes[k2]["lower"] <= v2) - assert (v2 <= self.G.nodes[k2]["upper"]) - # Check arrival times - for k1, v1 in self.prob.arrival_time.items(): - for k2, v2 in v1.items(): - assert (self.G.nodes[k2]["lower"] <= v2) - assert (v2 <= self.G.nodes[k2]["upper"]) + self.prob.check_arrival_time() + self.prob.check_departure_time() diff --git a/benchmarks/tests/test_cvrptw_solomon_range.py b/benchmarks/tests/test_cvrptw_solomon_range.py index 43c986e..f1bd9d3 100644 --- a/benchmarks/tests/test_cvrptw_solomon_range.py +++ b/benchmarks/tests/test_cvrptw_solomon_range.py @@ -1,10 +1,12 @@ from pytest import fixture +from time import time +import csv from vrpy import VehicleRoutingProblem from benchmarks.solomon_dataset import SolomonDataSet -params = list(range(7, 10)) +params = list(range(7, 70)) @fixture( @@ -16,17 +18,55 @@ def n(request): return request.param +REPS_LP = 1 +REPS_CSPY = 10 + + +def write_avg(n, times_cspy, iter_cspy, times_lp, iter_lp, name="cspy102fwdearly"): + def _avg(l): + return sum(l) / len(l) + + with open(f"benchmarks/results/{name}.csv", "a", newline="") as f: + writer_object = csv.writer(f) + writer_object.writerow( + [n, _avg(times_cspy), _avg(iter_cspy), _avg(times_lp), _avg(iter_lp)] + ) + f.close() + + class TestsSolomon: def test_subproblem(self, n): - data = SolomonDataSet(path="benchmarks/data/cvrptw/", - instance_name="C101.txt", - n_vertices=n) + data = SolomonDataSet( + path="benchmarks/data/cvrptw/", instance_name="C101.txt", n_vertices=n + ) self.G = data.G - self.prob = VehicleRoutingProblem(self.G, - load_capacity=data.max_load, - time_windows=True) - self.prob.solve(cspy=False) - best_value_lp = self.prob.best_value - self.prob.solve(cspy=True) - best_value_cspy = self.prob.best_value - assert int(best_value_lp) == int(best_value_cspy) + best_values_lp = None + lp_iter = [] + times_lp = [] + for r in range(REPS_LP): + prob = VehicleRoutingProblem( + self.G, load_capacity=data.max_load, time_windows=True + ) + start = time() + prob.solve(cspy=False) + best_value_lp = prob.best_value + times_lp.append(time() - start) + lp_iter.append(prob._iteration) + del prob + best_values_cspy = [] + times_cspy = [] + iter_cspy = [] + for r in range(REPS_CSPY): + prob = VehicleRoutingProblem( + self.G, load_capacity=data.max_load, time_windows=True + ) + start = time() + prob.solve(cspy=True, pricing_strategy="Exact") + times_cspy.append(time() - start) + best_values_cspy.append(prob.best_value) + iter_cspy.append(prob._iteration) + prob.check_arrival_time() + prob.check_departure_time() + del prob + assert all(best_value_lp == val_cspy for val_cspy in best_values_cspy) + write_avg(n, times_cspy, iter_cspy, times_lp, lp_iter) diff --git a/benchmarks/tests/test_examples.py b/benchmarks/tests/test_examples.py index 1bd249c..dec16a6 100644 --- a/benchmarks/tests/test_examples.py +++ b/benchmarks/tests/test_examples.py @@ -42,35 +42,39 @@ def setup(self): # Define VRP self.prob = VehicleRoutingProblem(self.G) - def test_cvrp_dive(self): + def test_cvrp_dive_lp(self): self.prob.load_capacity = 15 self.prob.solve(cspy=False, pricing_strategy="BestEdges1", dive=True) - sol_lp = self.prob.best_value + assert int(self.prob.best_value) == 6208 + + def test_cvrp_dive_cspy(self): + self.prob.load_capacity = 15 self.prob.solve(pricing_strategy="BestEdges1", dive=True) - sol_cspy = self.prob.best_value - assert int(sol_lp) == 6208 - assert int(sol_cspy) == 6208 + assert int(self.prob.best_value) == 6208 - def test_vrptw_dive(self): + def test_vrptw_dive_lp(self): self.prob.time_windows = True self.prob.solve(cspy=False, dive=True) - sol_lp = self.prob.best_value - self.prob.solve(dive=True) - sol_cspy = self.prob.best_value - assert int(sol_lp) == 6528 - assert int(sol_cspy) == 6528 + assert int(self.prob.best_value) == 6528 + + def test_vrptw_dive_cspy(self): + self.prob.time_windows = True + self.prob.solve(cspy=True, dive=True) + assert int(self.prob.best_value) == 6528 - def test_cvrpsdc_dive(self): + def test_cvrpsdc_dive_lp(self): self.prob.load_capacity = 15 self.prob.distribution_collection = True self.prob.solve(cspy=False, pricing_strategy="BestEdges1", dive=True) - sol_lp = self.prob.best_value + assert int(self.prob.best_value) == 6208 + + def test_cvrpsdc_dive_cspy(self): + self.prob.load_capacity = 15 + self.prob.distribution_collection = True self.prob.solve(pricing_strategy="BestEdges1", dive=True) - sol_cspy = self.prob.best_value - assert int(sol_lp) == 6208 - assert int(sol_cspy) == 6208 + assert int(self.prob.best_value) == 6208 - def test_pdp_dive(self): + def test_pdp_dive_lp(self): # Set demands and requests for (u, v) in PICKUPS_DELIVERIES: self.G.nodes[u]["request"] = v @@ -83,35 +87,39 @@ def test_pdp_dive(self): sol_lp = self.prob.best_value assert int(sol_lp) == 5980 - def test_cvrp(self): + def test_cvrp_lp(self): self.prob.load_capacity = 15 self.prob.solve(cspy=False, pricing_strategy="BestEdges1") - sol_lp = self.prob.best_value + assert int(self.prob.best_value) == 6208 + + def test_cvrp_cspy(self): + self.prob.load_capacity = 15 self.prob.solve(pricing_strategy="BestEdges1") - sol_cspy = self.prob.best_value - assert int(sol_lp) == 6208 - assert int(sol_cspy) == 6208 + assert int(self.prob.best_value) == 6208 - def test_vrptw(self): + def test_vrptw_lp(self): self.prob.time_windows = True self.prob.solve(cspy=False) - sol_lp = self.prob.best_value + assert int(self.prob.best_value) == 6528 + + def test_vrptw_cspy(self): + self.prob.time_windows = True self.prob.solve() - sol_cspy = self.prob.best_value - assert int(sol_lp) == 6528 - assert int(sol_cspy) == 6528 + assert int(self.prob.best_value) == 6528 - def test_cvrpsdc(self): + def test_cvrpsdc_lp(self): self.prob.load_capacity = 15 self.prob.distribution_collection = True self.prob.solve(cspy=False, pricing_strategy="BestEdges1") - sol_lp = self.prob.best_value + assert int(self.prob.best_value) == 6208 + + def test_cvrpsdc_cspy(self): + self.prob.load_capacity = 15 + self.prob.distribution_collection = True self.prob.solve(pricing_strategy="BestEdges1") - sol_cspy = self.prob.best_value - assert int(sol_lp) == 6208 - assert int(sol_cspy) == 6208 + assert int(self.prob.best_value) == 6208 - def test_pdp(self): + def test_pdp_lp(self): # Set demands and requests for (u, v) in PICKUPS_DELIVERIES: self.G.nodes[u]["request"] = v diff --git a/benchmarks/tests/test_issue101.py b/benchmarks/tests/test_issue101.py new file mode 100644 index 0000000..7a27ebf --- /dev/null +++ b/benchmarks/tests/test_issue101.py @@ -0,0 +1,19 @@ +from networkx import DiGraph, read_gpickle +from vrpy import VehicleRoutingProblem + + +class TestIssue101_large: + def setup(self): + G = read_gpickle("benchmarks/tests/graph_issue101") + self.prob = VehicleRoutingProblem(G, load_capacity=80) + self.prob.time_windows = True + + # def test_lp(self): + # self.prob.solve(cspy=False, solver="gurobi") + # self.prob.check_arrival_time() + # self.prob.check_departure_time() + + def test_cspy(self): + self.prob.solve(pricing_strategy="Exact") + self.prob.check_arrival_time() + self.prob.check_departure_time() diff --git a/examples/pdp.py b/examples/pdp.py index ef44f2e..6591557 100644 --- a/examples/pdp.py +++ b/examples/pdp.py @@ -21,8 +21,18 @@ if __name__ == "__main__": prob = VehicleRoutingProblem(G, load_capacity=6, pickup_delivery=True, num_stops=6) - prob.solve(cspy=False) + prob.solve(cspy=False, pricing_strategy="Exact") print(prob.best_value) print(prob.best_routes) + for (u, v) in PICKUPS_DELIVERIES: + found = False + for route in prob.best_routes.values(): + if u in route and v in route: + found = True + break + if not found: + print((u, v), "Not present") + assert False + print(prob.node_load) assert prob.best_value == 5980 diff --git a/tests/pytest.ini b/tests/pytest.ini index 5e0b9f0..a18e67e 100644 --- a/tests/pytest.ini +++ b/tests/pytest.ini @@ -1,5 +1,10 @@ [pytest] +log_cli = 1 +log_cli_level = INFO +log_cli_format = %(asctime)s [%(levelname)8s] %(message)s (%(filename)s:%(lineno)s) +log_cli_date_format=%Y-%m-%d %H:%M:%S + python_files = test_*.py #testpaths = tests/ filterwarnings = - ignore::DeprecationWarning \ No newline at end of file + ignore::DeprecationWarning diff --git a/tests/test_consistency.py b/tests/test_consistency.py index 763e17a..0bbdc08 100644 --- a/tests/test_consistency.py +++ b/tests/test_consistency.py @@ -53,9 +53,6 @@ def test_consistency_parameters(): G = DiGraph() G.add_edge("Source", "Sink", cost=1) prob = VehicleRoutingProblem(G, pickup_delivery=True) - # pickup delivery requires cspy=False - with pytest.raises(NotImplementedError): - prob.solve() # pickup delivery expects at least one request with pytest.raises(KeyError): prob.solve(cspy=False, pricing_strategy="Exact") @@ -88,9 +85,10 @@ def test_mixed_fleet_consistency(): prob.solve() G.edges["Source", "Sink"]["cost"] = [1, 2] with pytest.raises(ValueError): - prob = VehicleRoutingProblem( - G, mixed_fleet=True, load_capacity=[2, 4], fixed_cost=[4] - ) + prob = VehicleRoutingProblem(G, + mixed_fleet=True, + load_capacity=[2, 4], + fixed_cost=[4]) prob.solve() diff --git a/tests/test_issue101.py b/tests/test_issue101.py new file mode 100644 index 0000000..eaf3d14 --- /dev/null +++ b/tests/test_issue101.py @@ -0,0 +1,31 @@ +from networkx import DiGraph +from vrpy import VehicleRoutingProblem + + +class TestIssue101_small: + def setup(self): + self.G = DiGraph() + self.G.add_edge("Source", 1, cost=5) + self.G.add_edge("Source", 2, cost=5) + self.G.add_edge(1, "Sink", cost=5) + self.G.add_edge(2, "Sink", cost=5) + self.G.add_edge(1, 2, cost=1) + self.G.nodes[1]["lower"] = 0 + self.G.nodes[1]["upper"] = 20 + self.G.nodes[2]["lower"] = 0 + self.G.nodes[2]["upper"] = 20 + self.G.nodes[1]["service_time"] = 5 + self.G.nodes[2]["service_time"] = 5 + self.G.nodes[1]["demand"] = 8 + self.G.nodes[2]["demand"] = 8 + self.prob = VehicleRoutingProblem(self.G, load_capacity=10, time_windows=True) + + def test_cspy(self): + self.prob.solve() + self.prob.check_arrival_time() + self.prob.check_departure_time() + + def test_lp(self): + self.prob.solve(cspy=False) + self.prob.check_arrival_time() + self.prob.check_departure_time() diff --git a/tests/test_issue79.py b/tests/test_issue79.py index ad661af..bcbbd58 100644 --- a/tests/test_issue79.py +++ b/tests/test_issue79.py @@ -3,6 +3,7 @@ class TestIssue79: + def setup(self): G = DiGraph() G.add_edge("Source", 8, cost=0) @@ -24,13 +25,20 @@ def setup(self): G.nodes[6]["collect"] = 1 G.nodes[2]["collect"] = 1 G.nodes[5]["collect"] = 2 - self.prob = VehicleRoutingProblem( - G, load_capacity=15, distribution_collection=True - ) + self.prob = VehicleRoutingProblem(G, + load_capacity=15, + distribution_collection=True) - def test_node_load(self): + def test_node_load_cspy(self): self.prob.solve() assert self.prob.node_load[1][8] == 8 assert self.prob.node_load[1][6] == 5 assert self.prob.node_load[1][2] == 5 assert self.prob.node_load[1][5] == 5 + + def test_node_load_lp(self): + self.prob.solve(cspy=False) + assert self.prob.node_load[1][8] == 8 + assert self.prob.node_load[1][6] == 5 + assert self.prob.node_load[1][2] == 5 + assert self.prob.node_load[1][5] == 5 diff --git a/tests/test_issue99.py b/tests/test_issue99.py index 3736134..0ee8745 100644 --- a/tests/test_issue99.py +++ b/tests/test_issue99.py @@ -3,7 +3,7 @@ from vrpy import VehicleRoutingProblem -class TestIssue79: +class TestIssue99: def setup(self): distance_ = [ @@ -247,6 +247,12 @@ def setup(self): # Define VRP self.prob = VehicleRoutingProblem(G_, load_capacity=100) - def test_node_load(self): - self.prob.solve() + def test_lp(self): + self.prob.solve(cspy=False) + print(self.prob.best_routes) + assert self.prob.best_value == 829 + + def test_cspy(self): + self.prob.solve(cspy=True) + print(self.prob.best_routes) assert self.prob.best_value == 829 diff --git a/tests/test_toy.py b/tests/test_toy.py index 8e23278..5d7da96 100644 --- a/tests/test_toy.py +++ b/tests/test_toy.py @@ -41,7 +41,7 @@ def test_cspy_stops(self): ["Source", 4, 5, "Sink"], ] assert set(prob.best_routes_cost.values()) == {30, 40} - prob.solve(exact=False) + prob.solve() assert prob.best_value == 70 def test_cspy_stops_capacity(self): @@ -58,7 +58,7 @@ def test_cspy_stops_capacity_duration(self): with stop, capacity and duration constraints """ prob = VehicleRoutingProblem(self.G, num_stops=3, load_capacity=10, duration=62) - prob.solve(exact=False) + prob.solve() assert prob.best_value == 85 assert set(prob.best_routes_duration.values()) == {41, 62} assert prob.node_load[1]["Sink"] in [5, 10] @@ -85,16 +85,10 @@ def test_cspy_schedule(self): time_windows=True, ) prob.solve() - # Check departure times - for k1, v1 in prob.departure_time.items(): - for k2, v2 in v1.items(): - assert self.G.nodes[k2]["lower"] <= v2 - assert v2 <= self.G.nodes[k2]["upper"] - # Check arrival times - for k1, v1 in prob.arrival_time.items(): - for k2, v2 in v1.items(): - assert self.G.nodes[k2]["lower"] <= v2 - assert v2 <= self.G.nodes[k2]["upper"] + assert prob.departure_time[1]["Source"] == 0 + assert prob.arrival_time[1]["Sink"] in [41, 62] + prob.check_arrival_time() + prob.check_departure_time() ############### # subsolve lp # @@ -143,16 +137,8 @@ def test_LP_schedule(self): time_windows=True, ) prob.solve(cspy=False) - # Check departure times - for k1, v1 in prob.departure_time.items(): - for k2, v2 in v1.items(): - assert self.G.nodes[k2]["lower"] <= v2 - assert v2 <= self.G.nodes[k2]["upper"] - # Check arrival times - for k1, v1 in prob.arrival_time.items(): - for k2, v2 in v1.items(): - assert self.G.nodes[k2]["lower"] <= v2 - assert v2 <= self.G.nodes[k2]["upper"] + prob.check_arrival_time() + prob.check_departure_time() def test_LP_stops_elementarity(self): """Tests column generation procedure on toy graph""" @@ -246,7 +232,7 @@ def test_extend_preassignment(self): prob.solve(preassignments=routes) assert prob.best_value == 70 - def test_pick_up_delivery(self): + def test_pick_up_delivery_lp(self): self.G.nodes[2]["request"] = 5 self.G.nodes[2]["demand"] = 10 self.G.nodes[3]["demand"] = 10 @@ -263,6 +249,23 @@ def test_pick_up_delivery(self): prob.solve(pricing_strategy="Exact", cspy=False) assert prob.best_value == 65 + # def test_pick_up_delivery_cspy(self): + # self.G.nodes[2]["request"] = 5 + # self.G.nodes[2]["demand"] = 10 + # self.G.nodes[3]["demand"] = 10 + # self.G.nodes[3]["request"] = 4 + # self.G.nodes[4]["demand"] = -10 + # self.G.nodes[5]["demand"] = -10 + # self.G.add_edge(2, 5, cost=10) + # self.G.remove_node(1) + # prob = VehicleRoutingProblem( + # self.G, + # load_capacity=15, + # pickup_delivery=True, + # ) + # prob.solve(pricing_strategy="Exact", cspy=True) + # assert prob.best_value == 65 + def test_distribution_collection(self): self.G.nodes[1]["collect"] = 12 self.G.nodes[4]["collect"] = 1 diff --git a/vrpy/hyper_heuristic.py b/vrpy/hyper_heuristic.py index 012073d..7350c49 100644 --- a/vrpy/hyper_heuristic.py +++ b/vrpy/hyper_heuristic.py @@ -125,17 +125,14 @@ def pick_heuristic(self): pass # choose according to MAB maxval = max(self.heuristic_points.values()) - best_heuristics = [ - i for i, j in self.heuristic_points.items() if j == maxval - ] + best_heuristics = [i for i, j in self.heuristic_points.items() if j == maxval] if len(best_heuristics) == 1: self.current_heuristic = best_heuristics[0] else: self.current_heuristic = self.random_state.choice(best_heuristics) return self.current_heuristic - def update_scaling_factor(self, no_improvement_count: int, - no_improvement_iteration: int): + def update_scaling_factor(self, no_improvement_count: int): """ Implements Drake et al. (2012) @@ -180,8 +177,8 @@ def current_performance( elif self.performance_measure_type == "weighted_average": self._current_performance_wgtavr( - active_columns=active_columns, - new_objective_value=new_objective_value) + active_columns=active_columns, new_objective_value=new_objective_value + ) else: raise ValueError("performence_measure not set correctly!") self.set_current_objective(objective=new_objective_value) @@ -220,13 +217,11 @@ def reward(self, y, stagnated=True): x *= 0.9 return x - def update_parameters(self, iteration: int, no_improvement_count: int, - no_improvement_iteration: int): + def update_parameters(self, iteration: int, no_improvement_count: int): "Updates the high-level parameters" # measure time and add to weighted average self.iteration = iteration - self.update_scaling_factor(no_improvement_count, - no_improvement_iteration) + self.update_scaling_factor(no_improvement_count) # compute average of runtimes if self.performance_measure_type == "relative_improvement": self._update_params_relimp() @@ -245,40 +240,48 @@ def _compute_last_runtime(self): def _current_performance_relimp(self, produced_column: bool = False): # time measure if self.iteration > self.start_computing_average + 1: - self.d = max((self.average_runtime - self.last_runtime) / - self.average_runtime * 100, 0) - logger.debug("Resolve count %s, improvement %s", produced_column, - self.d) + self.d = max( + (self.average_runtime - self.last_runtime) / self.average_runtime * 100, + 0, + ) + logger.debug("Resolve count %s, improvement %s", produced_column, self.d) if self.d > self.d_max: self.d_max = self.d logger.debug( "Column produced, average runtime %s and last runtime %s", - self.average_runtime, self.last_runtime) + self.average_runtime, + self.last_runtime, + ) else: self.d = 0 - def _current_performance_wgtavr(self, - new_objective_value: float = None, - active_columns: dict = None): + def _current_performance_wgtavr( + self, new_objective_value: float = None, active_columns: dict = None + ): self.active_columns = active_columns # insert new runtime into sorted list bisect.insort(self.runtime_dist, self.last_runtime) self.objective_decrease[self.current_heuristic] += max( - self.current_objective_value - new_objective_value, 0) + self.current_objective_value - new_objective_value, 0 + ) self.total_objective_decrease += max( - self.current_objective_value - new_objective_value, 0) + self.current_objective_value - new_objective_value, 0 + ) # update quality values for heuristic in self.heuristic_options: if self.iterations[heuristic] > 0: self._update_exp(heuristic) - index = bisect.bisect(self.runtime_dist, - self.last_runtime_dict[heuristic]) - self.norm_runtime[heuristic] = (len(self.runtime_dist) - - index) / len(self.runtime_dist) + index = bisect.bisect( + self.runtime_dist, self.last_runtime_dict[heuristic] + ) + self.norm_runtime[heuristic] = (len(self.runtime_dist) - index) / len( + self.runtime_dist + ) if self.total_objective_decrease > 0: self.norm_objective_decrease[heuristic] = ( - self.objective_decrease[heuristic] / - self.total_objective_decrease) + self.objective_decrease[heuristic] + / self.total_objective_decrease + ) def _update_params_relimp(self): "Updates params for relative improvements performance measure" @@ -288,9 +291,10 @@ def _update_params_relimp(self): if reduced_n == 0: self.average_runtime = self.last_runtime else: - self.average_runtime = (self.average_runtime * (reduced_n - 1) / - (reduced_n) + 1 / - (reduced_n) * self.last_runtime) + self.average_runtime = ( + self.average_runtime * (reduced_n - 1) / (reduced_n) + + 1 / (reduced_n) * self.last_runtime + ) heuristic = self.current_heuristic # store old values old_q = self.q[heuristic] @@ -298,19 +302,22 @@ def _update_params_relimp(self): stagnated = old_q == 0 and old_n > 3 # average of improvements self.r[heuristic] = self.r[heuristic] * old_n / (old_n + 1) + 1 / ( - old_n + 1) * self.reward(self.d, stagnated=stagnated) + old_n + 1 + ) * self.reward(self.d, stagnated=stagnated) self.q[heuristic] = (old_q + self.r[heuristic]) / (old_n + 1) # compute heuristic points MAB-style for heuristic in self.heuristic_options: if self.iterations[heuristic] != 0: self._update_exp(heuristic) self.heuristic_points[heuristic] = ( - self.theta * self.q[heuristic] + - self.scaling_factor * self.exp[heuristic]) + self.theta * self.q[heuristic] + + self.scaling_factor * self.exp[heuristic] + ) def _update_exp(self, heuristic): - self.exp[heuristic] = sqrt(2 * log(sum(self.iterations.values())) / - self.iterations[heuristic]) + self.exp[heuristic] = sqrt( + 2 * log(sum(self.iterations.values())) / self.iterations[heuristic] + ) def _update_params_wgtavr(self): "Updates params for Weighted average performance measure" @@ -326,10 +333,11 @@ def _update_params_wgtavr(self): if sum_exp != 0: norm_spread = self.exp[heuristic] / sum_exp self.q[heuristic] = ( - self.weight_col_basic * active_i / active + - self.weight_runtime * norm_runtime + - self.weight_obj * self.norm_objective_decrease[heuristic] + - self.weight_col_total * active_i / total_added) + self.weight_col_basic * active_i / active + + self.weight_runtime * norm_runtime + + self.weight_obj * self.norm_objective_decrease[heuristic] + + self.weight_col_total * active_i / total_added + ) self.heuristic_points[heuristic] = self.theta * self.q[ - heuristic] + self.weight_spread * norm_spread * (1 - - self.theta) + heuristic + ] + self.weight_spread * norm_spread * (1 - self.theta) diff --git a/vrpy/master_solve_pulp.py b/vrpy/master_solve_pulp.py index c9988c3..3bcbe97 100644 --- a/vrpy/master_solve_pulp.py +++ b/vrpy/master_solve_pulp.py @@ -7,6 +7,7 @@ from vrpy.masterproblem import _MasterProblemBase from vrpy.restricted_master_heuristics import _DivingHeuristic + logger = logging.getLogger(__name__) @@ -50,11 +51,13 @@ def solve(self, relax, time_limit): if pulp.LpStatus[self.prob.status] != "Optimal": raise Exception("problem " + str(pulp.LpStatus[self.prob.status])) - if relax: - for r in self.routes: - val = pulp.value(self.y[r.graph["name"]]) - if val > 0.1: - logger.debug("route %s selected %s" % (r.graph["name"], val)) + # This logging takes time + # if relax: + # for r in self.routes: + # val = pulp.value(self.y[r.graph["name"]]) + # if val > 0.1: + # logger.debug("route %s selected %s" % + # (r.graph["name"], val)) duals = self.get_duals() logger.debug("duals : %s" % duals) return duals, self.prob.objective.value() @@ -123,16 +126,13 @@ def get_total_cost_and_routes(self, relax: bool): for r in self.routes: val = pulp.value(self.y[r.graph["name"]]) if val is not None and val > 0: - logger.debug( - "%s cost %s load %s" - % ( - shortest_path(r, "Source", "Sink"), - r.graph["cost"], - sum(self.G.nodes[v]["demand"] for v in r.nodes()), - ) - ) - best_routes.append(r) + # This logging takes time + # logger.debug("%s cost %s load %s" % ( + # shortest_path(r, "Source", "Sink"), + # r.graph["cost"], + # sum(self.G.nodes[v]["demand"] for v in r.nodes()), + # )) if self.drop_penalty: self.dropped_nodes = [ v for v in self.drop if pulp.value(self.drop[v]) > 0.5 @@ -178,7 +178,9 @@ def _solve(self, relax: bool, time_limit: Optional[int]): and "depot_to" not in self.G.nodes[node] ): for const in self.prob.constraints: - # Modify the self.prob object (the self.set_covering_constrs object cannot be modified (?)) + # Modify the self.prob object (the + # self.set_covering_constrs object cannot be modified + # (?)) if "visit_node" in const: self.prob.constraints[const].sense = pulp.LpConstraintEQ if ( @@ -189,6 +191,10 @@ def _solve(self, relax: bool, time_limit: Optional[int]): self.dummy[node].cat = pulp.LpInteger # Set route variables to integer for var in self.y.values(): + # Disallow routes that visit multiple nodes + if "non" in var.name: + var.upBound = 0 + var.lowBound = 0 var.cat = pulp.LpInteger # Force vehicle bound artificial variable to 0 for var in self.dummy_bound.values(): @@ -303,17 +309,20 @@ def _add_set_covering_constraints(self): ) def _add_route_selection_variable(self, route): + # Added path with the raw path as using `nx.add_path` then + # `graph.nodes()` + # removes repeated nodes, so the coefficients are not correct. + if "path" in route.graph: + nodes = [n for n in route.graph["path"] if n not in ["Source", "Sink"]] + else: + nodes = [n for n in route.nodes() if n not in ["Source", "Sink"]] self.y[route.graph["name"]] = pulp.LpVariable( "y{}".format(route.graph["name"]), lowBound=0, upBound=1, cat=pulp.LpContinuous, e=( - pulp.lpSum( - self.set_covering_constrs[r] - for r in route.nodes() - if r not in ["Source", "Sink"] - ) + pulp.lpSum(self.set_covering_constrs[n] for n in nodes) + pulp.lpSum( self.vehicle_bound_constrs[k] for k in range(len(self.num_vehicles)) diff --git a/vrpy/subproblem.py b/vrpy/subproblem.py index a25b294..308609e 100644 --- a/vrpy/subproblem.py +++ b/vrpy/subproblem.py @@ -198,7 +198,7 @@ def remove_edges_3(self, beta): 4. Remove all edges that do not belong to these paths """ # Normalize weights - max_weight = max([self.G.edges[i, j]["weight"] for (i, j) in self.G.edges()]) + max_weight = max(self.G.edges[i, j]["weight"] for (i, j) in self.G.edges()) min_weight = min(self.G.edges[i, j]["weight"] for (i, j) in self.G.edges()) for edge in self.G.edges(data=True): edge[2]["pos_weight"] = ( diff --git a/vrpy/subproblem_cspy.py b/vrpy/subproblem_cspy.py index 31c1ef7..789f044 100644 --- a/vrpy/subproblem_cspy.py +++ b/vrpy/subproblem_cspy.py @@ -1,8 +1,10 @@ import logging from math import floor + from numpy import zeros from networkx import DiGraph, add_path from cspy import BiDirectional, REFCallback + from vrpy.subproblem import _SubProblemBase logger = logging.getLogger(__name__) @@ -14,8 +16,14 @@ class _MyREFCallback(REFCallback): Based on Righini and Salani (2006). """ - def __init__(self, max_res, time_windows, distribution_collection, T, - resources): + def __init__( + self, + max_res, + time_windows, + distribution_collection, + T, + resources, + ): REFCallback.__init__(self) # Set attributes for use in REF functions self._max_res = max_res @@ -28,42 +36,29 @@ def __init__(self, max_res, time_windows, distribution_collection, T, self._source_id = None self._sink_id = None - def REF_fwd(self, cumul_res, tail, head, edge_res, partial_path, - cumul_cost): + def REF_fwd(self, cumul_res, tail, head, edge_res, partial_path, cumul_cost): new_res = list(cumul_res) i, j = tail, head # stops / monotone resource new_res[0] += 1 # load new_res[1] += self._sub_G.nodes[j]["demand"] + # time # Service times theta_i = self._sub_G.nodes[i]["service_time"] - # theta_j = self._sub_G.nodes[j]["service_time"] - # theta_t = self._sub_G.nodes[self._sink_id]["service_time"] # Travel times travel_time_ij = self._sub_G.edges[i, j]["time"] - # try: - # travel_time_jt = self._sub_G.edges[j, self._sink_id]["time"] - # except KeyError: - # travel_time_jt = 0 # Time windows # Lower a_j = self._sub_G.nodes[j]["lower"] - # a_t = self._sub_G.nodes[self._sink_id]["lower"] # Upper b_j = self._sub_G.nodes[j]["upper"] - # b_t = self._sub_G.nodes[self._sink_id]["upper"] new_res[2] = max(new_res[2] + theta_i + travel_time_ij, a_j) # time-window feasibility resource if not self._time_windows or (new_res[2] <= b_j): - # and new_res[2] < self._T - a_j - theta_j and - # a_t <= new_res[2] + travel_time_jt + theta_t <= b_t): - # if not self._time_windows or ( - # new_res[2] <= b_j and new_res[2] < self._T - a_j - theta_j and - # a_t <= new_res[2] + travel_time_jt + theta_t <= b_t): new_res[3] = 0 else: new_res[3] = 1 @@ -72,12 +67,11 @@ def REF_fwd(self, cumul_res, tail, head, edge_res, partial_path, # Pickup new_res[4] += self._sub_G.nodes[j]["collect"] # Delivery - new_res[5] = max(new_res[5] + self._sub_G.nodes[j]["demand"], - new_res[4]) + new_res[5] = max(new_res[5] + self._sub_G.nodes[j]["demand"], new_res[4]) + return new_res - def REF_bwd(self, cumul_res, tail, head, edge_res, partial_path, - cumul_cost): + def REF_bwd(self, cumul_res, tail, head, edge_res, partial_path, cumul_cost): new_res = list(cumul_res) i, j = tail, head @@ -85,37 +79,23 @@ def REF_bwd(self, cumul_res, tail, head, edge_res, partial_path, new_res[0] -= 1 # load new_res[1] += self._sub_G.nodes[i]["demand"] + # Get relevant service times (thetas) and travel time # Service times theta_i = self._sub_G.nodes[i]["service_time"] theta_j = self._sub_G.nodes[j]["service_time"] - # theta_s = self._sub_G.nodes[self._source_id]["service_time"] # Travel times travel_time_ij = self._sub_G.edges[i, j]["time"] - # try: - # travel_time_si = self._sub_G.edges[self._source_id, i]["time"] - # except KeyError: - # travel_time_si = 0 # Lower time windows - # a_i = self._sub_G.nodes[i]["lower"] a_j = self._sub_G.nodes[j]["lower"] - # a_s = self._sub_G.nodes[self._source_id]["lower"] # Upper time windows b_i = self._sub_G.nodes[i]["upper"] b_j = self._sub_G.nodes[j]["upper"] - # b_s = self._sub_G.nodes[self._source_id]["upper"] - new_res[2] = max(new_res[2] + theta_j + travel_time_ij, - self._T - b_i - theta_i) + new_res[2] = max(new_res[2] + theta_j + travel_time_ij, self._T - b_i - theta_i) # time-window feasibility if not self._time_windows or (new_res[2] <= self._T - a_j - theta_j): - # and new_res[2] < self._T - a_i - theta_i and - # a_s <= new_res[2] + theta_s + travel_time_si <= b_s): - # if not self._time_windows or ( - # new_res[2] <= self._T - a_j and - # new_res[2] < self._T - a_i - theta_i and - # a_s <= new_res[2] + theta_s + travel_time_si <= b_s): new_res[3] = 0 else: new_res[3] = 1 @@ -123,9 +103,8 @@ def REF_bwd(self, cumul_res, tail, head, edge_res, partial_path, if self._distribution_collection: # Delivery new_res[5] += new_res[5] + self._sub_G.nodes[i]["demand"] - # Pickup - new_res[4] = max(new_res[5], - new_res[4] + self._sub_G.nodes[i]["collect"]) + # Pick up + new_res[4] = max(new_res[5], new_res[4] + self._sub_G.nodes[i]["collect"]) return new_res def REF_join(self, fwd_resources, bwd_resources, tail, head, edge_res): @@ -172,11 +151,11 @@ class _SubProblemCSPY(_SubProblemBase): Inherits problem parameters from `SubproblemBase` """ - def __init__(self, *args, exact): + def __init__(self, *args, elementary): """Initializes resources.""" # Pass arguments to base super(_SubProblemCSPY, self).__init__(*args) - self.exact = exact + self.elementary = elementary # Resource names self.resources = [ "stops/mono", @@ -191,13 +170,13 @@ def __init__(self, *args, exact): # Default lower and upper bounds self.min_res = [0] * len(self.resources) # Add upper bounds for mono, stops, load and time, and time windows - total_demand = sum( - [self.sub_G.nodes[v]["demand"] for v in self.sub_G.nodes()]) + total_demand = sum([self.sub_G.nodes[v]["demand"] for v in self.sub_G.nodes()]) self.max_res = [ floor(len(self.sub_G.nodes()) / 2), # stop/mono total_demand, # load - sum([self.sub_G.edges[u, v]["time"] for u, v in self.sub_G.edges() - ]), # time + sum( + [self.sub_G.edges[u, v]["time"] for u, v in self.sub_G.edges()] + ), # time 1, # time windows total_demand, # pickup total_demand, # deliver @@ -208,6 +187,10 @@ def __init__(self, *args, exact): # Initialize max feasible arrival time self.T = 0 self.total_cost = None + # Average length of a path + self._avg_path_len = 1 + # Iteration counter + self._iters = 1 # @profile def solve(self, time_limit): @@ -215,10 +198,10 @@ def solve(self, time_limit): Solves the subproblem with cspy. Time limit is reduced by 0.5 seconds as a safety window. - Resolves until: - 1. heuristic algorithm gives a new route (column with -ve reduced cost); - 2. exact algorithm gives a new route; - 3. neither heuristic nor exact give a new route. + Resolves at most twice: + 1. using elementary = False, + 2. using elementary = True, and threshold, if a route has already been + found previously. """ if not self.run_subsolve: return self.routes, False @@ -231,28 +214,46 @@ def solve(self, time_limit): more_routes = False my_callback = self.get_REF() - while True: - if self.exact: - alg = BiDirectional( - self.sub_G, - self.max_res, - self.min_res, - direction="both", - time_limit=time_limit - 0.5 if time_limit else None, - elementary=True, - REF_callback=my_callback, + direction = ( + "forward" + if ( + self.time_windows + or self.pickup_delivery + or self.distribution_collection + ) + else "both" + ) + # Run only twice: Once with `elementary=False` check if route already + # exists. + + s = ( + [False, True] + if (not self.distribution_collection and not self.elementary) + else [True] + ) + for elementary in s: + if elementary: + # Use threshold if non-elementary (safe-guard against large + # instances) + thr = self._avg_path_len * min( + self.G.edges[i, j]["weight"] for (i, j) in self.G.edges() ) else: - alg = BiDirectional( - self.sub_G, - self.max_res, - self.min_res, - threshold=-1, - direction="both", - time_limit=time_limit - 0.5 if time_limit else None, - elementary=True, - REF_callback=my_callback, - ) + thr = None + logger.debug( + f"Solving subproblem using elementary={elementary}, threshold={thr}, direction={direction}" + ) + alg = BiDirectional( + self.sub_G, + self.max_res, + self.min_res, + threshold=thr, + direction=direction, + time_limit=time_limit - 0.5 if time_limit else None, + elementary=elementary, + REF_callback=my_callback, + # pickup_delivery_pairs=self.pickup_delivery_pairs, + ) # Pass processed graph if my_callback is not None: @@ -263,19 +264,27 @@ def solve(self, time_limit): logger.debug("subproblem") logger.debug("cost = %s", alg.total_cost) logger.debug("resources = %s", alg.consumed_resources) + if alg.total_cost is not None and alg.total_cost < -(1e-3): - more_routes = True - self.add_new_route(alg.path) - logger.debug("new route %s", alg.path) - logger.debug("reduced cost = %s", alg.total_cost) - logger.debug("real cost = %s", self.total_cost) - break - # If not already solved exactly - elif not self.exact: - # Solve exactly from here on - self.exact = True - # Solved heuristically and exactly and no more routes - # Or time out + new_route = self.create_new_route(alg.path) + logger.debug(alg.path) + path_len = len(alg.path) + if not any( + list(new_route.edges()) == list(r.edges()) for r in self.routes + ): + more_routes = True + self.routes.append(new_route) + self.total_cost = new_route.graph["cost"] + logger.debug("reduced cost = %s", alg.total_cost) + logger.debug("real cost = %s", self.total_cost) + if path_len > 2: + self._avg_path_len += ( + path_len - self._avg_path_len + ) / self._iters + self._iters += 1 + break + else: + logger.info("Route already found, finding elementary one") else: break return self.routes, more_routes @@ -298,33 +307,34 @@ def formulate(self): # Time windows feasibility self.max_res[3] = 0 # Maximum feasible arrival time - # for v in self.sub_G.nodes(): - # print("node = ", v, "lb = ", self.sub_G.nodes[v]["lower"], - # "ub = ", self.sub_G.nodes[v]["upper"]) - - self.T = max(self.sub_G.nodes[v]["upper"] + - self.sub_G.nodes[v]["service_time"] + - self.sub_G.edges[v, "Sink"]["time"] - for v in self.sub_G.predecessors("Sink")) + self.T = max( + self.sub_G.nodes[v]["upper"] + + self.sub_G.nodes[v]["service_time"] + + self.sub_G.edges[v, "Sink"]["time"] + for v in self.sub_G.predecessors("Sink") + ) if self.load_capacity and self.distribution_collection: self.max_res[4] = self.load_capacity[self.vehicle_type] self.max_res[5] = self.load_capacity[self.vehicle_type] + if self.pickup_delivery: + self.max_res[6] = 0 - def add_new_route(self, path): + def create_new_route(self, path): """Create new route as DiGraph and add to pool of columns""" - route_id = len(self.routes) + 1 - new_route = DiGraph(name=route_id) + e = "elem" if len(set(path)) == len(path) else "non-elem" + route_id = "{}_{}".format(len(self.routes) + 1, e) + new_route = DiGraph(name=route_id, path=path) add_path(new_route, path) - self.total_cost = 0 + total_cost = 0 for (i, j) in new_route.edges(): edge_cost = self.sub_G.edges[i, j]["cost"][self.vehicle_type] - self.total_cost += edge_cost + total_cost += edge_cost new_route.edges[i, j]["cost"] = edge_cost if i != "Source": self.routes_with_node[i].append(new_route) - new_route.graph["cost"] = self.total_cost + new_route.graph["cost"] = total_cost new_route.graph["vehicle_type"] = self.vehicle_type - self.routes.append(new_route) + return new_route def add_max_stops(self): """Updates maximum number of stops.""" @@ -357,7 +367,7 @@ def add_max_duration(self): self.sub_G.edges[i, j]["res_cost"][2] = travel_time def get_REF(self): - if self.time_windows or self.distribution_collection: + if self.time_windows or self.distribution_collection or self.pickup_delivery: # Use custom REF return _MyREFCallback( self.max_res, diff --git a/vrpy/subproblem_lp.py b/vrpy/subproblem_lp.py index 84f60c8..73a38a3 100644 --- a/vrpy/subproblem_lp.py +++ b/vrpy/subproblem_lp.py @@ -41,7 +41,7 @@ def solve(self, time_limit, exact=True): logger.debug("Objective %s" % pulp.value(self.prob.objective)) if ( pulp.value(self.prob.objective) is not None - and pulp.value(self.prob.objective) < -(10 ** -3) + and pulp.value(self.prob.objective) < -(10**-3) and pulp.LpStatus[self.prob.status] not in ["Infeasible"] ) or (exact == False and pulp.LpStatus[self.prob.status] in ["Optimal", ""]): more_routes = True diff --git a/vrpy/vrp.py b/vrpy/vrp.py index 2d27e39..11a0690 100644 --- a/vrpy/vrp.py +++ b/vrpy/vrp.py @@ -1,14 +1,19 @@ +import sys import logging from time import time from typing import List, Union + from networkx import DiGraph, shortest_path # draw_networkx + from vrpy.greedy import _Greedy -from vrpy.master_solve_pulp import _MasterSolvePulp +from vrpy.schedule import _Schedule from vrpy.subproblem_lp import _SubProblemLP from vrpy.subproblem_cspy import _SubProblemCSPY +from vrpy.hyper_heuristic import _HyperHeuristic +from vrpy.master_solve_pulp import _MasterSolvePulp from vrpy.subproblem_greedy import _SubProblemGreedy +from vrpy.preprocessing import get_num_stops_upper_bound from vrpy.clarke_wright import _ClarkeWright, _RoundTrip -from vrpy.schedule import _Schedule from vrpy.checks import ( check_arguments, check_consistency, @@ -20,12 +25,10 @@ check_clarke_wright_compatibility, check_preassignments, ) -from vrpy.preprocessing import get_num_stops_upper_bound -from vrpy.hyper_heuristic import _HyperHeuristic logger = logging.getLogger(__name__) -logging.basicConfig(level=logging.INFO) +logging.basicConfig(level=logging.INFO, stream=sys.stdout) class VehicleRoutingProblem: @@ -120,8 +123,8 @@ def __init__( self._solver: str = None self._time_limit: int = None self._pricing_strategy: str = None - self._exact: bool = None self._cspy: bool = None + self._elementary: bool = None self._dive: bool = None self._greedy: bool = None self._max_iter: int = None @@ -153,7 +156,7 @@ def solve( preassignments=None, pricing_strategy="BestEdges1", cspy=True, - exact=False, + elementary=False, time_limit=None, solver="cbc", dive=False, @@ -185,10 +188,12 @@ def solve( cspy (bool, optional): True if cspy is used for subproblem. Defaults to True. - exact (bool, optional): - True if only cspy's exact algorithm is used to generate columns. - Otherwise, heuristics will be used until they produce +ve - reduced cost columns, after which the exact algorithm is used. + elementary (bool, optional): + True if only cspy's elementary algorithm is to be used. + Otherwise, a mix is used: only use elementary when a route is + repeated. In this case, also a threshold is used to avoid + running for too long. If dive=True then this is also forced to + be True. Defaults to False. time_limit (int, optional): Maximum number of seconds allowed for solving (for finding columns). @@ -224,8 +229,8 @@ def solve( self._solver = solver self._time_limit = time_limit self._pricing_strategy = pricing_strategy - self._exact = exact self._cspy = cspy + self._elementary = elementary if not dive else True self._dive = False self._greedy = greedy self._max_iter = max_iter @@ -233,7 +238,6 @@ def solve( self._heuristic_only = heuristic_only if self._pricing_strategy == "Hyper": self.hyper_heuristic = _HyperHeuristic() - self._start_time = time() if preassignments: self._preassignments = preassignments @@ -367,17 +371,21 @@ def arrival_time(self): for j in range(1, len(route)): tail = route[j - 1] head = route[j] - arrival[i][head] = min( - max( - arrival[i][tail] - + self._H.nodes[tail]["service_time"] - + self._H.edges[tail, head]["time"], - self._H.nodes[head]["lower"], - ), - self._H.nodes[head]["upper"], + arrival[i][head] = max( + arrival[i][tail] + + self._H.nodes[tail]["service_time"] + + self._H.edges[tail, head]["time"], + self._H.nodes[head]["lower"], ) return arrival + def check_arrival_time(self): + # Check arrival times + for k1, v1 in self.arrival_time.items(): + for k2, v2 in v1.items(): + assert self.G.nodes[k2]["lower"] <= v2 + assert v2 <= self.G.nodes[k2]["upper"] + @property def departure_time(self): """ @@ -395,12 +403,20 @@ def departure_time(self): for j in range(1, len(route)): tail = route[j - 1] head = route[j] - departure[i][head] = min( - arrival[i][head] + self._H.nodes[head]["service_time"], - self._H.nodes[head]["upper"], + departure[i][head] = ( + arrival[i][head] + self._H.nodes[head]["service_time"] ) return departure + def check_departure_time(self): + # Check departure times + for k1, v1 in self.departure_time.items(): + for k2, v2 in v1.items(): + assert self.G.nodes[k2]["lower"] <= v2 + # Upper TW should not be checked as lower + service_time > upper + # in many cases. Also not enforced at the subproblem level. + # assert v2 <= self.G.nodes[k2]["upper"] + @property def schedule(self): """If Periodic CVRP, returns a dict with keys a day number and values @@ -508,6 +524,7 @@ def _solve(self, dive, solver): self._best_value, self._best_routes_as_graphs, ) = self.masterproblem.get_total_cost_and_routes(relax=False) + self._post_process(solver) def _column_generation(self): @@ -543,7 +560,9 @@ def _find_columns(self): logger.info("iteration %s, %.6s" % (self._iteration, relaxed_cost)) pricing_strategy = self._get_next_pricing_strategy(relaxed_cost) + # TODO: parallel # One subproblem per vehicle type + for vehicle in range(self._vehicle_types): # Solve pricing problem with randomised greedy algorithm if ( @@ -593,7 +612,6 @@ def _find_columns(self): self._no_improvement += 1 else: self._no_improvement = 0 - self._no_improvement_iteration = self._iteration if not self._dive: self._lower_bound.append(relaxed_cost) @@ -674,7 +692,7 @@ def _attempt_solve_best_edges1(self, vehicle=None, duals=None, route=None): alpha, ) self.routes, self._more_routes = subproblem.solve( - self._get_time_remaining(), + self._get_time_remaining() ) more_columns = self._more_routes if more_columns: @@ -695,7 +713,6 @@ def _attempt_solve_best_edges2(self, vehicle=None, duals=None, route=None): ) self.routes, self._more_routes = subproblem.solve( self._get_time_remaining(), - # exact=False, ) more_columns = self._more_routes if more_columns: @@ -717,7 +734,7 @@ def _get_next_pricing_strategy(self, relaxed_cost): self._pricing_strategy == "Hyper" and self._no_improvement != self._run_exact ): - self._no_improvement_iteration = self._iteration + self._no_improvement = self._iteration if self._iteration == 0: pricing_strategy = "BestPaths" self.hyper_heuristic.init(relaxed_cost) @@ -726,7 +743,7 @@ def _get_next_pricing_strategy(self, relaxed_cost): self._update_hyper_heuristic(relaxed_cost) pricing_strategy = self.hyper_heuristic.pick_heuristic() elif self._no_improvement == self._run_exact: - self._no_improvement = 0 + # self._no_improvement = 0 pricing_strategy = "Exact" else: pricing_strategy = self._pricing_strategy @@ -740,9 +757,7 @@ def _update_hyper_heuristic(self, relaxed_cost: float): active_columns=best_paths_freq, ) self.hyper_heuristic.move_acceptance() - self.hyper_heuristic.update_parameters( - self._iteration, self._no_improvement, self._no_improvement_iteration - ) + self.hyper_heuristic.update_parameters(self._iteration, self._no_improvement) def _get_time_remaining(self, mip: bool = False): """ @@ -808,7 +823,7 @@ def _def_subproblem( self.distribution_collection, pricing_strategy, pricing_parameter, - exact=self._exact, + elementary=self._elementary, ) else: # As LP @@ -1051,7 +1066,7 @@ def _prune_graph(self): self._remove_infeasible_arcs_time_windows() def _set_zero_attributes(self): - """ Sets attr = 0 if missing """ + """Sets attr = 0 if missing""" for v in self.G.nodes(): for attribute in [ @@ -1076,14 +1091,14 @@ def _set_zero_attributes(self): self.G.nodes[v][attribute] = 1 def _set_time_to_zero_if_missing(self): - """ Sets time = 0 if missing """ + """Sets time = 0 if missing""" for (i, j) in self.G.edges(): for attribute in ["time"]: if attribute not in self.G.edges[i, j]: self.G.edges[i, j][attribute] = 0 def _readjust_sink_time_windows(self): - """ Readjusts Sink time windows """ + """Readjusts Sink time windows""" if self.G.nodes["Sink"]["upper"] == 0: self.G.nodes["Sink"]["upper"] = max(