google / or-tools

Google's Operations Research tools:
https://developers.google.com/optimization/
Apache License 2.0
11.27k stars 2.13k forks source link

SlackVar for multi-depot example does not work as expected #1127

Closed quinlivanb closed 3 years ago

quinlivanb commented 5 years ago

The example below is a very simple vrp with capacity constraint, requiring one return journey to the depot. The truck should be able to travel to nodes 2 and 3 and accumulate a load of 20/25. It should then drop this at the depot, which has a demand of -25 and a SlackVar range of 0 25 before returning to node 4 and then back to the actual depot.

Example of a possible solution:

0 Load(0) -> 2 Load(10) -> 3 Load(20) -> 1 Load(0) -> 4 Load(20) -> 0 Load(20).

However the solver does not find this (or any) solution. The SlackVar at the dummy depot should allow the above route but it is not working expected.

Also note, if you set the capacity to 30 then it will solve the problem as it can now return to the dummy depot with the exact load required.

from __future__ import print_function
import math
from ortools.constraint_solver import pywrapcp
from ortools.constraint_solver import routing_enums_pb2

def euclid_distance(x1, y1, x2, y2):
    # Euclidean distance between points.
    dist = math.sqrt((x1 - x2)**2 + (y1 - y2)**2)
    return dist

def create_distance_matrix(locations):
    # Create the distance matrix.
    size = len(locations)
    dist_matrix = {}

    for from_node in range(size):
        dist_matrix[from_node] = {}
        for to_node in range(size):
            x1 = locations[from_node][0]
            y1 = locations[from_node][1]
            x2 = locations[to_node][0]
            y2 = locations[to_node][1]
            dist_matrix[from_node][to_node] = euclid_distance(x1, y1, x2, y2)
    return dist_matrix

def create_data_model():
    """Creates the data for the example."""
    data = {}
    # Array of distances between locations.
    _locations = [
        (0, 0),
        (0, 0),
        (-50, 50),
        (0, 100),
        (50, 50)]
    capacities = 25
    _distances = create_distance_matrix(_locations)
    demands = [
        0,
        -capacities,
        10,
        10,
        20]
    data["distances"] = _distances
    data["num_locations"] = len(_distances)
    data["num_vehicles"] = 1
    data["depot"] = 0
    data["demands"] = demands
    data["vehicle_capacities"] = capacities
    return data

def create_distance_callback(data):
    """Creates callback to return distance between points."""
    distances = data["distances"]

    def distance_callback(from_node, to_node):
        """Returns the manhattan distance between the two nodes"""
        return distances[from_node][to_node]
    return distance_callback

def create_demand_callback(data):
    """Creates callback to get demands at each location."""
    def demand_callback(from_node, to_node):
        del to_node
        return data["demands"][from_node]
    return demand_callback

def add_capacity_constraints(routing, data, demand_callback):
    """Adds capacity constraint"""
    capacity = "Capacity"
    routing.AddDimension(
        demand_callback,
        0,  # null capacity slack
        data["vehicle_capacities"],  # vehicle maximum capacities
        True,  # start cumul to zero
        capacity)
    # Add Slack for reseting to zero unload depot nodes.
    # e.g. vehicle with load 10/15 arrives at node 1 (depot unload)
    # so we have CumulVar = 10(current load) + -15(unload) + 5(slack) = 0.
    capacity_dimension = routing.GetDimensionOrDie(capacity)
    for node_index in [1]:
        index = routing.NodeToIndex(node_index)
        capacity_dimension.SlackVar(index).SetRange(0, data["vehicle_capacities"])

def print_solution(data, routing, assignment):
    """Print routes on console."""
    total_dist = 0
    for vehicle_id in range(data["num_vehicles"]):
        index = routing.Start(vehicle_id)
        plan_output = 'Route for vehicle {0}:\n'.format(vehicle_id)
        route_dist = 0
        route_load = 0
        while not routing.IsEnd(index):
            node_index = routing.IndexToNode(index)
            next_node_index = routing.IndexToNode(assignment.Value(routing.NextVar(index)))
            route_dist += routing.GetArcCostForVehicle(node_index, next_node_index, vehicle_id)
            route_load += data["demands"][node_index]
            plan_output += ' {0} Load({1}) -> '.format(node_index, route_load)
            index = assignment.Value(routing.NextVar(index))

        node_index = routing.IndexToNode(index)
        total_dist += route_dist
        plan_output += ' {0} Load({1})\n'.format(node_index, route_load)
        plan_output += 'Distance of the route: {0}m\n'.format(route_dist)
        plan_output += 'Load of the route: {0}\n'.format(route_load)
        print(plan_output)
    print('Total Distance of all routes: {0}m'.format(total_dist))

def main():
    """Entry point of the program"""
    # Instantiate the data problem.
    data = create_data_model()
    # Create Routing Model
    routing = pywrapcp.RoutingModel(
        data["num_locations"],
        data["num_vehicles"],
        data["depot"])
    # Define weight of each edge
    distance_callback = create_distance_callback(data)
    routing.SetArcCostEvaluatorOfAllVehicles(distance_callback)
    # Add Capacity constraint
    demand_callback = create_demand_callback(data)
    add_capacity_constraints(routing, data, demand_callback)
    # Setting first solution heuristic (cheapest addition).
    search_parameters = pywrapcp.RoutingModel.DefaultSearchParameters()
    search_parameters.first_solution_strategy = (
        routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC)
    # Solve the problem.
    assignment = routing.SolveWithParameters(search_parameters)
    if assignment:
        print_solution(data, routing, assignment)
    else:
        print('No solution has been found')

if __name__ == '__main__':
    main()
lkillora commented 5 years ago

SlackVar does not seem to work at all. Below is a simple example with a depot and one node. They have demands of 0 and -5 respectively. Now set the SlackVar to any range, and it won't allow the solution 0Load(0) -> 1Load(-5) -> 0Load(-5). Of course, if you set the capacity slack from null (0) to 5, then it will work, but for this, no SlackVar is needed. In summary, SlackVar does not seem to affect the outcomes in any way.

    from __future__ import print_function
    import math
    from ortools.constraint_solver import pywrapcp
    from ortools.constraint_solver import routing_enums_pb2

    def euclid_distance(x1, y1, x2, y2):
        # Euclidean distance between points.
        dist = math.sqrt((x1 - x2)**2 + (y1 - y2)**2)
        return dist

    def create_distance_matrix(locations):
        # Create the distance matrix.
        size = len(locations)
        dist_matrix = {}

        for from_node in range(size):
            dist_matrix[from_node] = {}
            for to_node in range(size):
                x1 = locations[from_node][0]
                y1 = locations[from_node][1]
                x2 = locations[to_node][0]
                y2 = locations[to_node][1]
                dist_matrix[from_node][to_node] = euclid_distance(x1, y1, x2, y2)
        return dist_matrix

    def create_data_model():
        """Creates the data for the example."""
        data = {}
        # Array of distances between locations.
        _locations = [
            (0, 0),
            (0, 0)]
        capacities = 25
        _distances = create_distance_matrix(_locations)
        demands = [
            0,
            -5]
        data["distances"] = _distances
        data["num_locations"] = len(_distances)
        data["num_vehicles"] = 1
        data["depot"] = 0
        data["demands"] = demands
        data["vehicle_capacities"] = capacities
        return data

    def create_distance_callback(data):
        """Creates callback to return distance between points."""
        distances = data["distances"]

        def distance_callback(from_node, to_node):
            """Returns the manhattan distance between the two nodes"""
            return distances[from_node][to_node]
        return distance_callback

    def create_demand_callback(data):
        """Creates callback to get demands at each location."""
        def demand_callback(from_node, to_node):
            del to_node
            return data["demands"][from_node]
        return demand_callback

    def add_capacity_constraints(routing, data, demand_callback):
        """Adds capacity constraint"""
        capacity = "Capacity"
        routing.AddDimension(
            demand_callback,
            0,  # null capacity slack
            data["vehicle_capacities"],  # vehicle maximum capacities
            True,  # start cumul to zero
            capacity)
        # Add Slack
        # e.g. vehicle with load 0/5 arrives at node 1 (depot unload)
        # so we have CumulVar = 0(current load) + -5(unload) + 5(slack) = 0.
        capacity_dimension = routing.GetDimensionOrDie(capacity)
        for node_index in [0, 1]:
            index = routing.NodeToIndex(node_index)
            capacity_dimension.SlackVar(index).SetRange(0, 5)

    def print_solution(data, routing, assignment):
        """Print routes on console."""
        total_dist = 0
        for vehicle_id in range(data["num_vehicles"]):
            index = routing.Start(vehicle_id)
            plan_output = 'Route for vehicle {0}:\n'.format(vehicle_id)
            route_dist = 0
            route_load = 0
            while not routing.IsEnd(index):
                node_index = routing.IndexToNode(index)
                next_node_index = routing.IndexToNode(assignment.Value(routing.NextVar(index)))
                route_dist += routing.GetArcCostForVehicle(node_index, next_node_index, vehicle_id)
                route_load += data["demands"][node_index]
                plan_output += ' {0} Load({1}) -> '.format(node_index, route_load)
                index = assignment.Value(routing.NextVar(index))

            node_index = routing.IndexToNode(index)
            total_dist += route_dist
            plan_output += ' {0} Load({1})\n'.format(node_index, route_load)
            plan_output += 'Distance of the route: {0}m\n'.format(route_dist)
            plan_output += 'Load of the route: {0}\n'.format(route_load)
            print(plan_output)
        print('Total Distance of all routes: {0}m'.format(total_dist))

    def main():
        """Entry point of the program"""
        # Instantiate the data problem.
        data = create_data_model()
        # Create Routing Model
        routing = pywrapcp.RoutingModel(
            data["num_locations"],
            data["num_vehicles"],
            data["depot"])
        # Define weight of each edge
        distance_callback = create_distance_callback(data)
        routing.SetArcCostEvaluatorOfAllVehicles(distance_callback)
        # Add Capacity constraint
        demand_callback = create_demand_callback(data)
        add_capacity_constraints(routing, data, demand_callback)
        # Setting first solution heuristic (cheapest addition).
        search_parameters = pywrapcp.RoutingModel.DefaultSearchParameters()
        search_parameters.first_solution_strategy = (
            routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC)
        # Solve the problem.
        assignment = routing.SolveWithParameters(search_parameters)
        if assignment:
            print_solution(data, routing, assignment)
        else:
            print('No solution has been found')

    if __name__ == '__main__':
        main()
Mizux commented 5 years ago

CumulVar must be positive or zero, if you want to use negative value you need to have positive SlackVar to compensate...

Mizux commented 5 years ago

did you try the sample: https://github.com/google/or-tools/blob/master/ortools/constraint_solver/samples/cvrp_reload.py

see #862

quinlivanb commented 5 years ago

Above slackVar is set to the range of 0 and capacity (25). The vrp_reload example has the same issue. A truck will only return the the dummy depots if it's capacity is exactly full. If you use the above locations, demands and capacity in the vrp_reload example it will also fail to find a route.

quinlivanb commented 5 years ago

Ikillora... Mizux is correct, cumulvar cannot be negative in your example so the expected solution with Slack would be

0Load(0) -> 1Load(0) -> 0Load(0).

But the fact that this isn't working also shows that this slackVar method shown in vrp_reload does not work correctly. Mizux, would really appreciate your help on this issue.

lkillora commented 5 years ago

@Mizux, what work is SlackVar doing in the routing.cc file and, specifically, in the FeasibleRoute function (provided below)?

    bool FeasibleRoute(const std::vector<int>& route, int64 route_cumul,
                     int dimension_index) {
    const RoutingDimension& dimension = *dimensions_[dimension_index];
    std::vector<int>::const_iterator it = route.begin();
    int64 cumul = route_cumul;
    while (it != route.end()) {
      const int previous = *it;
      const int64 cumul_previous = cumul;
      gtl::InsertOrDie(&(new_possible_cumuls_[dimension_index]), previous,
                       cumul_previous);
      ++it;
      if (it == route.end()) {
        return true;
      }
      const int next = *it;
      int64 available_from_previous =
          cumul_previous + dimension.GetTransitValue(previous, next, 0);
      int64 available_cumul_next =
          std::max(cumuls_[dimension_index][next], available_from_previous);

      const int64 slack = available_cumul_next - available_from_previous;
      if (slack > dimension.SlackVar(previous)->Max()) {
        available_cumul_next =
            available_from_previous + dimension.SlackVar(previous)->Max();
      }

      if (available_cumul_next > dimension.CumulVar(next)->Max()) {
        return false;
      }
      if (available_cumul_next <= cumuls_[dimension_index][next]) {
        return true;
      }
      cumul = available_cumul_next;
    }
    return true;
    }
quinlivanb commented 5 years ago

@Mizux Hey, do you have any more input on this? We are still blocked by this issue. I've tested on every available example using capacity slack and they all have the same issue. The truck will only return to dummy depot when the capacity of the node matches the load in the truck exactly.

in the cvrp_reload.py example, removing the following line has no effect. And If you increase the capacity to 17 it will start dropping nodes.

capacity_dimension.SlackVar(index).SetRange(0, vehicle_capacity)

I hope you can help with this, thanks!

270ajay commented 5 years ago

@quinlivanb

there is another way to model that: instead of setting 0 slack in the dimension and providing a range of slack for dummy node, we can set max slack in the dimension, and set 0 slack for non-dummy nodes.

if you update the add_capacity_constraints to this, it works:

def add_capacity_constraints(routing, data, demand_callback):
    """Adds capacity constraint"""
    capacity = "Capacity"
    routing.AddDimension(
        demand_callback,
        data["vehicle_capacities"],  # null capacity slack
        data["vehicle_capacities"],  # vehicle maximum capacities
        True,  # start cumul to zero
        capacity)
    # Add Slack for reseting to zero unload depot nodes.
    # e.g. vehicle with load 10/15 arrives at node 1 (depot unload)
    # so we have CumulVar = 10(current load) + -15(unload) + 5(slack) = 0.
    capacity_dimension = routing.GetDimensionOrDie(capacity)
    for node_index in range(routing.nodes()):
        index = routing.NodeToIndex(node_index)
        if node_index not in [1]:
            capacity_dimension.SlackVar(index).SetValue(0)
        routing.AddVariableMinimizedByFinalizer(capacity_dimension.CumulVar(index))

you may also need to update print_solution function Output:

Route for vehicle 0:
 0 Load(0) ->  4 Load(20) ->  1 Load(-5) ->  2 Load(5) ->  3 Load(15) ->  0 Load(15)
Distance of the route: 380m
Load of the route: 15

Total Distance of all routes: 380m

Downside of this way of modeling: Sometimes you may see vehicle visiting dummy node/reloading when it may not be required. However, it is not very frequent.

jiangmouren commented 5 years ago

Guys, any updates on this?

jmarca commented 5 years ago

Poking around, perhaps this line is the problem: https://github.com/google/or-tools/blob/811ab90009d4886e6cd7c75730ccf250c190e4e1/ortools/constraint_solver/routing.cc#L5797-L5802

So in the original way, slack_max is zero in the definition of the dimension. So there are two issues. First, the slack is never taken into account in the transit_expr. And second, according to expressions.cc, line 2576 or so, calling "SetRange" on a IntConst variable doesn't quite do what its name might imply:

class IntConst : public IntVar {
 public:
  IntConst(Solver* const s, int64 value, const std::string& name = "")
      : IntVar(s, name), value_(value) {}
  ~IntConst() override {}
  ...
  void SetRange(int64 l, int64 u) override {
    if (l > value_ || u < value_) {
      solver()->Fail();
    }
  }

In the rearranging way suggested by @270ajay the slackmax is not zero, so the slacks[i] array are all IntVar's with a range from 0 to max, which ends up being a DomainIntVar, not IntConst, and the transit expression actually cares about the value of the slack variable. The implementation of SetRange for DomainIntVar actually does do what you think it does.

So I think this is a bug. Not sure how to fix it though, other than to change the example to be like @270ajay suggestion?

jmarca commented 5 years ago

Note that you can't just say

if (slack_max == 0) {
    slacks_[i] = solver->MakeIntVar(0, slack_max, slack_name);
    transit_expr = solver->MakeSum(slacks_[i], transit_expr);
}

because MakeIntVar right away has the optimization that any range of 0 is immediately an IntConst value:

IntVar* Solver::MakeIntVar(int64 min, int64 max, const std::string& name) {
  if (min == max) {
    return MakeIntConst(min, name);
  }
  if (min == 0 && max == 1) {
    return RegisterIntVar(RevAlloc(new ConcreteBooleanVar(this, name)));
  } else if (CapSub(max, min) == 1) {
    const std::string inner_name = "inner_" + name;
    return RegisterIntVar(
        MakeSum(RevAlloc(new ConcreteBooleanVar(this, inner_name)), min)
            ->VarWithName(name));
  } else {
    return RegisterIntVar(RevAlloc(new DomainIntVar(this, min, max, name)));
  }
}

Which again means you can't go back in after the fact and change the min and max of this constant value.

amih77 commented 5 years ago

Hello guys, any news about this?

I have a slightly different use-case where I need recharge when the vehicle cumulative distance exceeds some threshold (so the capacities are defined over the arcs rather than the nodes).

It's described here https://github.com/google/or-tools/issues/1624

hherbertb commented 4 years ago

Hi, are there any updates on this issue?

vc-void commented 4 years ago

@quinlivanb

there is another way to model that: instead of setting 0 slack in the dimension and providing a range of slack for dummy node, we can set max slack in the dimension, and set 0 slack for non-dummy nodes.

if you update the add_capacity_constraints to this, it works:

def add_capacity_constraints(routing, data, demand_callback):
    """Adds capacity constraint"""
    capacity = "Capacity"
    routing.AddDimension(
        demand_callback,
        data["vehicle_capacities"],  # null capacity slack
        data["vehicle_capacities"],  # vehicle maximum capacities
        True,  # start cumul to zero
        capacity)
    # Add Slack for reseting to zero unload depot nodes.
    # e.g. vehicle with load 10/15 arrives at node 1 (depot unload)
    # so we have CumulVar = 10(current load) + -15(unload) + 5(slack) = 0.
    capacity_dimension = routing.GetDimensionOrDie(capacity)
    for node_index in range(routing.nodes()):
        index = routing.NodeToIndex(node_index)
        if node_index not in [1]:
            capacity_dimension.SlackVar(index).SetValue(0)
        routing.AddVariableMinimizedByFinalizer(capacity_dimension.CumulVar(index))

you may also need to update print_solution function Output:

Route for vehicle 0:
 0 Load(0) ->  4 Load(20) ->  1 Load(-5) ->  2 Load(5) ->  3 Load(15) ->  0 Load(15)
Distance of the route: 380m
Load of the route: 15

Total Distance of all routes: 380m

Downside of this way of modeling: Sometimes you may see vehicle visiting dummy node/reloading when it may not be required. However, it is not very frequent.

This workaround worked for me in java, although I had to increase the distance between depots to avoid the mentioned downside.

hherbertb commented 4 years ago

Thank you for your quick response. I also saw this one, but I hoped there is something new available.

Based on your experience, what is a good distance value between dummy nodes and the rest of the nodes?

vc-void commented 4 years ago

That depends really. I am still testing and instead of distances i have time windows. For now i have set it so that going from a depot to another depot will hit one of my time window constraints, invalidating routes with 2 depots in a row.

qesbalcha commented 4 years ago

Any update on this? For the time being, I am just adding some offset to the vehicle capacity compared to the depot capacity... And this works but, the amount of offset is just intuitive.

jmarca commented 4 years ago

Updating my previous comment, the link should now be https://github.com/google/or-tools/blob/master/ortools/constraint_solver/routing.cc#L5797

This just came up again in the response to issue #1891 and on the follow-up on the mailing list.

I think maybe the code should fail rather than succeeding when trying to call SetRange(a,b) on an IntConst.

At https://github.com/google/or-tools/blob/811ab90009d4886e6cd7c75730ccf250c190e4e1/ortools/constraint_solver/expressions.cc#L2571-L2575
, maybe the code should be

  void SetRange(int64 l, int64 u) override {
      solver()->Fail();
  }

Why? Because if you call set range on an IntConst, the caller probably thinks that the full range is available, but instead it will just be the single value.

Or promote the IntConst to IntVar

But silently succeeding is leading to mystery bugs