always-cautious/ortools-clone
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

graph: update dag constrained shortest path

Browse Source
This commit is contained in:
Corentin Le Molgat
2024-03-25 18:20:28 +01:00
parent f1f380d985
commit e666eda009
3 changed files with 783 additions and 251 deletions
Download Patch File Download Diff File Expand all files Collapse all files

4
ortools/graph/dag_constrained_shortest_path.cc
View File

@@ -78,8 +78,8 @@ PathWithLength ConstrainedShortestPathsOnDag(
std::vector<NodeIndex> destinations = {destination};
ConstrainedShortestPathsOnDagWrapper<util::StaticGraph<>>
constrained_shortest_path_on_dag(&graph, &arc_lengths, &arc_resources,
&(*topological_order), &sources,
&destinations, &max_resources);
*topological_order, sources,
destinations, &max_resources);
PathWithLength path_with_length =
constrained_shortest_path_on_dag.RunConstrainedShortestPathOnDag();

836
ortools/graph/dag_constrained_shortest_path.h
View File

@@ -14,6 +14,8 @@
#ifndef OR_TOOLS_GRAPH_DAG_CONSTRAINED_SHORTEST_PATH_H_
#define OR_TOOLS_GRAPH_DAG_CONSTRAINED_SHORTEST_PATH_H_
#include <stdbool.h>
#include <cmath>
#include <concepts>
#include <limits>
@@ -25,13 +27,25 @@
#include "absl/strings/str_format.h"
#include "absl/types/span.h"
#include "ortools/graph/dag_shortest_path.h"
#include "ortools/graph/graph.h"
namespace operations_research {
// This library provides APIs to compute the constrained shortest path (CSP) on
// a given directed acyclic graph (DAG) with resources on each arc. A CSP is a
// shortest path on a DAG which does not exceed a set of maximum resources
// consumption. The algorithm is exponential and has no guarantee to finish.
// consumption. The algorithm is exponential and has no guarantee to finish. It
// is based on bi-drectionnal search. First is a forward pass from the source to
// nodes “somewhere in the middle” to generate forward labels, just as the
// onedirectional labeling algorithm we discussed; then a symmetric backward
// pass from the destination generates backward labels; and finally at each node
// with both forward and backward labels, it joins any pair of labels to form a
// feasible complete path. Intuitively, the number of labels grows exponentially
// with the number of arcs in the path. The overall number of labels are then
// expected to be smaller with shorter paths. For DAG with a topological
// ordering, we can pick any node (usually right in the middle) as a *midpoint*
// to stop each pass at. Then labels can be joined at only one half of the nodes
// by considering all edges between each half.
//
// In the DAG, multiple arcs between the same pair of nodes is allowed. However,
// self-loop arcs are not allowed.
@@ -114,70 +128,156 @@ class ConstrainedShortestPathsOnDagWrapper {
//
// Validity of arcs and topological order are DCHECKed.
//
// If the number of labels in memory exceeds `max_num_created_labels` at any
// point in the algorithm, it returns the best path found so far, most
// particularly the empty path if none were found.
// If the number of labels in memory exceeds `max_num_created_labels / 2` at
// any point in each pass of the algorithm, new labels are not generated
// anymore and it returns the best path found so far, most particularly the
// empty path if none were found.
//
// SUBTLE: You can modify the graph, the arc lengths and resources, the
// topological order, sources, destinations or the maximum resource between
// calls to the `RunConstrainedShortestPathOnDag()` function. That's fine.
// Doing so will obviously invalidate the result API of the last constrained
// shortest path run, which could return an upper bound, junk, or crash.
// IMPORTANT: You cannot modify anything except `arc_lengths` between calls to
// the `RunConstrainedShortestPathOnDag()` function.
// TODO(b/330285675): Change API to copy everything internally except the
// lengths.
ConstrainedShortestPathsOnDagWrapper(
const GraphType* graph, const std::vector<double>* arc_lengths,
const std::vector<std::vector<double>>* arc_resources,
const std::vector<NodeIndex>* topological_order,
const std::vector<NodeIndex>* sources,
const std::vector<NodeIndex>* destinations,
absl::Span<const NodeIndex> topological_order,
absl::Span<const NodeIndex> sources,
absl::Span<const NodeIndex> destinations,
const std::vector<double>* max_resources,
int max_num_created_labels = 1e9);
// Returns {+inf, {}, {}} if there is no constrained path of finite length
// from one node in `sources` to one node in `destinations`.
// wihtin resources constraints from one node in `sources` to one node in
// `destinations`.
PathWithLength RunConstrainedShortestPathOnDag();
// For benchmarking and informational purposes, returns the number of labels
// generated in the call of `RunConstrainedShortestPathOnDag()`.
int label_count() const {
return lengths_from_sources_[FORWARD].size() +
lengths_from_sources_[BACKWARD].size();
}
private:
// Returns the list of all the arcs of the shortest path from one node in
// `sources` ending by the arc from a given `label_index` if and only if
// `label_index` is between 0 and `labels_from_sources_.size() - 1`.
std::vector<ArcIndex> BestArcPathEndingWith(int label_index) const;
enum Direction {
FORWARD = 0,
BACKWARD = 1,
};
inline static Direction Reverse(Direction d) {
return d == FORWARD ? BACKWARD : FORWARD;
}
// A LabelPair includes the `length` of a path that can be constructed by
// merging the paths from two *linkable* labels corresponding to
// `label_index`.
struct LabelPair {
double length = 0.0;
int label_index[2];
};
void RunHalfConstrainedShortestPathOnDag(
const GraphType& reverse_graph, absl::Span<const double> arc_lengths,
absl::Span<const std::vector<double>> arc_resources,
absl::Span<const std::vector<double>> min_arc_resources,
absl::Span<const double> max_resources, int max_num_created_labels,
std::vector<double>& lengths_from_sources,
std::vector<std::vector<double>>& resources_from_sources,
std::vector<ArcIndex>& incoming_arc_indices_from_sources,
std::vector<int>& first_label, std::vector<int>& num_labels);
// Returns the arc index linking two nodes from each pass forming the best
// path. Returns -1 if no better path than the one found from
// `best_label_pair` is found.
ArcIndex MergeHalfRuns(
const GraphType& graph, absl::Span<const double> arc_lengths,
absl::Span<const std::vector<double>> arc_resources,
absl::Span<const double> max_resources,
const std::vector<NodeIndex> sub_node_indices[2],
const std::vector<double> lengths_from_sources[2],
const std::vector<std::vector<double>> resources_from_sources[2],
const std::vector<int> first_label[2],
const std::vector<int> num_labels[2], LabelPair& best_label_pair);
// Returns the path as list of arc indices that starts from a node in
// `sources` (if `direction` iS FORWARD) or `destinations` (if `direction` is
// BACKWARD) and ends in node represented by `best_label_index`.
std::vector<ArcIndex> ArcPathTo(
int best_label_index, const GraphType& reverse_graph,
absl::Span<const double> arc_lengths,
absl::Span<const double> lengths_from_sources,
absl::Span<const ArcIndex> incoming_arc_indices_from_sources,
absl::Span<const int> first_label,
absl::Span<const int> num_labels) const;
// Returns the list of all the nodes implied by a given `arc_path`.
std::vector<NodeIndex> NodePathImpliedBy(
const std::vector<ArcIndex>& arc_path) const;
std::vector<NodeIndex> NodePathImpliedBy(absl::Span<const ArcIndex> arc_path,
const GraphType& graph) const;
static constexpr double kTolerance = 1e-6;
const GraphType* const graph_;
const std::vector<double>* const arc_lengths_;
const std::vector<std::vector<double>>* const arc_resources_;
const std::vector<NodeIndex>* const topological_order_;
const std::vector<NodeIndex>* const sources_;
const std::vector<NodeIndex>* const destinations_;
const std::vector<double>* const max_resources_;
int max_num_created_labels_;
absl::Span<const NodeIndex> sources_;
absl::Span<const NodeIndex> destinations_;
const int num_resources_;
std::vector<bool> node_is_source_;
std::vector<bool> node_is_destination_;
// Using the fact that the graph is a DAG, we can disregard any node that
// comes after the last destination (based on the topological order).
std::vector<bool> node_is_after_last_destination_;
// Data for reverse graph.
GraphType reverse_graph_;
std::vector<ArcIndex> reverse_inverse_arc_permutation_;
// Data about *reachable* sub-graphs split in two for bidirectional search.
// Reachable nodes are nodes that can be reached given the resources
// constraints, i.e., for each resource, the sum of the minimum resource to
// get to a node from a node in `sources` and to get from a node to a node in
// `destinations` should be less than the maximum resource. Reachable arcs are
// arcs linking reachable nodes.
//
// `sub_reverse_graph_[dir]` is the reachable sub-graph split in *half* with
// an additional linked to sources (resp. destinations) for the forward (resp.
// backward) direction. For the forward (resp. backward) direction, nodes are
// indexed using the original (resp. reverse) topological order.
GraphType sub_reverse_graph_[2];
std::vector<std::vector<double>> sub_arc_resources_[2];
// `sub_full_arc_indices_[dir]` has size `sub_reverse_graph_[dir].num_arcs()`
// such that `sub_full_arc_indices_[dir][sub_arc] = arc` where `sub_arc` is
// the arc in the reachable sub-graph for direction `dir` (i.e.
// `sub_reverse_graph[dir]`) and `arc` is the arc in the original graph (i.e.
// `graph`).
std::vector<NodeIndex> sub_full_arc_indices_[2];
// `sub_node_indices_[dir]` has size `graph->num_nodes()` such that
// `sub_node_indices[dir][node] = sub_node` where `node` is the node in the
// original graph (i.e. `graph`) and `sub_node` is the node in the reachable
// sub-graph for direction `dir` (i.e. `sub_reverse_graph[dir]`) and -1 if
// `node` is not present in reachable sub-graph.
std::vector<NodeIndex> sub_node_indices_[2];
// `sub_is_source_[dir][sub_dir]` has size
// `sub_reverse_graph_[dir].num_nodes()` such that
// `sub_is_source_[dir][sub_dir][sub_node]` is true if `sub_node` is a node in
// the reachable sub-graph for direction `dir` (i.e. `sub_reverse_graph[dir]`)
// which is a source (resp. destination) is `sub_dir` is FORWARD (resp.
// BACKWARD).
std::vector<bool> sub_is_source_[2][2];
// `sub_min_arc_resources_[dir]` has size `max_resources->size()` and
// `sub_min_arc_resources_[dir][r]`, `sub_reverse_graph_[dir].num_nodes()`
// such that `sub_min_arc_resources_[dir][r][sub_node]` is the minimum of
// resource r needed to get to a destination (resp. come from a source) if
// `dir` is FORWARD (resp. BACKWARD).
std::vector<std::vector<double>> sub_min_arc_resources_[2];
// Maximum number of labels created for each sub-graph.
int max_num_created_labels_[2];
// Data about the last call of the RunConstrainedShortestPathOnDag()
// function. A Label includes the cumulative length, resources and the
// previous arc used in the path to get to this node.
struct Label {
double length;
// TODO(b/315786885): Optimize resources in Label struct.
std::vector<double> resources;
ArcIndex incoming_arc;
};
// A label is present in `labels_from_sources_` if and only if it is feasible
// with respect to all resources.
std::vector<Label> labels_from_sources_;
std::vector<int> node_first_label_;
std::vector<int> node_num_labels_;
// function. A path is only added to the following vectors if and only if
// it is feasible with respect to all resources.
// A Label includes the cumulative length, resources and the previous arc used
// in the path to get to this node.
// Instead of having a single vector of `Label` objects (cl/590819865), we
// split them into 3 vectors of more fundamental types as this improves
// push_back operations and memory release.
std::vector<double> lengths_from_sources_[2];
std::vector<std::vector<double>> resources_from_sources_[2];
std::vector<ArcIndex> incoming_arc_indices_from_sources_[2];
std::vector<int> node_first_label_[2];
std::vector<int> node_num_labels_[2];
};
std::vector<int> GetInversePermutation(absl::Span<const int> permutation);
@@ -194,30 +294,31 @@ ConstrainedShortestPathsOnDagWrapper<GraphType>::
ConstrainedShortestPathsOnDagWrapper(
const GraphType* graph, const std::vector<double>* arc_lengths,
const std::vector<std::vector<double>>* arc_resources,
const std::vector<NodeIndex>* topological_order,
const std::vector<NodeIndex>* sources,
const std::vector<NodeIndex>* destinations,
absl::Span<const NodeIndex> topological_order,
absl::Span<const NodeIndex> sources,
absl::Span<const NodeIndex> destinations,
const std::vector<double>* max_resources, int max_num_created_labels)
: graph_(graph),
arc_lengths_(arc_lengths),
arc_resources_(arc_resources),
topological_order_(topological_order),
max_resources_(max_resources),
sources_(sources),
destinations_(destinations),
max_resources_(max_resources),
max_num_created_labels_(max_num_created_labels) {
num_resources_(max_resources->size()) {
CHECK(graph_ != nullptr);
CHECK(arc_lengths_ != nullptr);
CHECK(arc_resources_ != nullptr);
CHECK(topological_order_ != nullptr);
CHECK(sources_ != nullptr);
CHECK(destinations_ != nullptr);
CHECK(!sources_.empty());
CHECK(!destinations_.empty());
CHECK(max_resources_ != nullptr);
CHECK(!max_resources_->empty())
<< "max_resources cannot be empty. Use "
"ortools/graph/dag_shortest_path.h instead";
if (DEBUG_MODE) {
CHECK_EQ(arc_lengths_->size(), graph_->num_arcs());
CHECK_EQ(arc_resources_->size(), max_resources_->size());
for (const std::vector<double>& arcs_resource : *arc_resources_) {
CHECK_EQ(arcs_resource.size(), graph_->num_arcs());
CHECK_EQ(arc_lengths->size(), graph->num_arcs());
CHECK_EQ(arc_resources->size(), max_resources->size());
for (absl::Span<const double> arcs_resource : *arc_resources) {
CHECK_EQ(arcs_resource.size(), graph->num_arcs());
for (const double arc_resource : arcs_resource) {
CHECK(arc_resource >= 0 &&
arc_resource != std::numeric_limits<double>::infinity() &&
@@ -225,53 +326,192 @@ ConstrainedShortestPathsOnDagWrapper<GraphType>::
<< absl::StrFormat("resource cannot be negative nor +inf nor NaN");
}
}
for (const double arc_length : *arc_lengths_) {
for (const double arc_length : *arc_lengths) {
CHECK(arc_length != -std::numeric_limits<double>::infinity() &&
!std::isnan(arc_length))
<< absl::StrFormat("length cannot be -inf nor NaN");
}
CHECK_OK(TopologicalOrderIsValid(*graph_, *topological_order_))
CHECK_OK(TopologicalOrderIsValid(*graph, topological_order))
<< "Invalid topological order";
for (const double max_resource : *max_resources_) {
for (const double max_resource : *max_resources) {
CHECK(max_resource >= 0 &&
max_resource != std::numeric_limits<double>::infinity() &&
!std::isnan(max_resource))
<< absl::StrFormat(
"max_resource cannot be negative not +inf nor NaN");
}
}
node_is_source_.resize(graph_->num_nodes());
for (const NodeIndex source : *sources_) {
node_is_source_[source] = true;
}
node_is_destination_.resize(graph_->num_nodes());
for (const NodeIndex destination : *destinations_) {
node_is_destination_[destination] = true;
}
node_is_after_last_destination_.resize(graph_->num_nodes());
bool last_destination_found = false;
for (int i = graph_->num_nodes() - 1; i >= 0; --i) {
const NodeIndex node = (*topological_order_)[i];
if (node_is_destination_[node]) {
last_destination_found = true;
std::vector<bool> is_source(graph->num_nodes(), false);
for (const NodeIndex source : sources) {
is_source[source] = true;
}
for (const NodeIndex destination : destinations) {
CHECK(!is_source[destination])
<< "A node cannot be both a source and destination";
}
node_is_after_last_destination_[node] = !last_destination_found;
}
// Reverse graph.
reverse_graph_ = GraphType(graph_->num_nodes(), graph_->num_arcs());
for (ArcIndex arc_index = 0; arc_index < graph_->num_arcs(); ++arc_index) {
reverse_graph_.AddArc(graph_->Head(arc_index), graph_->Tail(arc_index));
// Full graphs.
const GraphType* full_graph[2];
const std::vector<std::vector<double>>* full_arc_resources[2];
absl::Span<const NodeIndex> full_topological_order[2];
absl::Span<const NodeIndex> full_sources[2];
// Forward.
const int num_nodes = graph->num_nodes();
const int num_arcs = graph->num_arcs();
full_graph[FORWARD] = graph;
full_arc_resources[FORWARD] = arc_resources;
full_topological_order[FORWARD] = topological_order;
full_sources[FORWARD] = sources;
// Backward.
GraphType full_backward_graph(num_nodes, num_arcs);
for (ArcIndex arc_index = 0; arc_index < num_arcs; ++arc_index) {
full_backward_graph.AddArc(graph->Head(arc_index), graph->Tail(arc_index));
}
std::vector<ArcIndex> permutation;
reverse_graph_.Build(&permutation);
reverse_inverse_arc_permutation_ = GetInversePermutation(permutation);
std::vector<ArcIndex> full_permutation;
full_backward_graph.Build(&full_permutation);
const std::vector<ArcIndex> full_inverse_arc_indices =
GetInversePermutation(full_permutation);
std::vector<std::vector<double>> backward_arc_resources(num_resources_);
for (int r = 0; r < num_resources_; ++r) {
backward_arc_resources[r] = (*arc_resources)[r];
util::Permute(full_permutation, &backward_arc_resources[r]);
}
std::vector<NodeIndex> full_backward_topological_order;
full_backward_topological_order.reserve(num_nodes);
for (int i = num_nodes - 1; i >= 0; --i) {
full_backward_topological_order.push_back(topological_order[i]);
}
full_graph[BACKWARD] = &full_backward_graph;
full_arc_resources[BACKWARD] = &backward_arc_resources;
full_topological_order[BACKWARD] = full_backward_topological_order;
full_sources[BACKWARD] = destinations;
// Memory allocation is done here and only once in order to avoid reallocation
// at each call of `RunConstrainedShortestPathOnDag()` for better performance.
node_first_label_.resize(graph_->num_nodes());
node_num_labels_.resize(graph_->num_nodes());
// Get the minimum resources sources -> node and node -> destination for each
// node.
std::vector<std::vector<double>> full_min_arc_resources[2];
for (const Direction dir : {FORWARD, BACKWARD}) {
full_min_arc_resources[dir].reserve(num_resources_);
std::vector<double> full_arc_resource = full_arc_resources[dir]->front();
ShortestPathsOnDagWrapper<GraphType> shortest_paths_on_dag(
full_graph[dir], &full_arc_resource, full_topological_order[dir]);
for (int r = 0; r < num_resources_; ++r) {
full_arc_resource = (*(full_arc_resources[dir]))[r];
shortest_paths_on_dag.RunShortestPathOnDag(full_sources[dir]);
full_min_arc_resources[dir].push_back(shortest_paths_on_dag.LengthTo());
}
}
// Get reachable subgraph.
std::vector<NodeIndex> sub_topological_order;
sub_topological_order.reserve(num_nodes);
for (const NodeIndex node_index : topological_order) {
bool is_reachable = true;
for (int r = 0; r < num_resources_; ++r) {
if (full_min_arc_resources[FORWARD][r][node_index] +
full_min_arc_resources[BACKWARD][r][node_index] >
(*max_resources)[r]) {
is_reachable = false;
break;
}
}
if (is_reachable) {
sub_topological_order.push_back(node_index);
}
}
const int reachable_node_count = sub_topological_order.size();
// Split sub-graphs and related information.
const int mid_index = reachable_node_count / 2;
// We split the number of labels evenly between each search (+1 for the
// additional source node)
// TODO(329206892): Could be optimized depending on the split and number of
// arcs in the sub-graph.
max_num_created_labels_[FORWARD] = max_num_created_labels / 2;
max_num_created_labels_[BACKWARD] =
max_num_created_labels - max_num_created_labels_[FORWARD];
for (const Direction dir : {FORWARD, BACKWARD}) {
absl::Span<const NodeIndex> const sub_nodes =
dir == FORWARD
? absl::MakeSpan(sub_topological_order).subspan(0, mid_index)
: absl::MakeSpan(sub_topological_order)
.subspan(mid_index, reachable_node_count - mid_index);
sub_node_indices_[dir].assign(num_nodes, -1);
sub_min_arc_resources_[dir].resize(num_resources_);
for (int r = 0; r < num_resources_; ++r) {
sub_min_arc_resources_[dir][r].resize(sub_nodes.size());
}
for (NodeIndex i = 0; i < sub_nodes.size(); ++i) {
const NodeIndex sub_node_index =
dir == FORWARD ? i : sub_nodes.size() - 1 - i;
sub_node_indices_[dir][sub_nodes[i]] = sub_node_index;
for (int r = 0; r < num_resources_; ++r) {
sub_min_arc_resources_[dir][r][sub_node_index] =
full_min_arc_resources[Reverse(dir)][r][sub_nodes[i]];
}
}
// IMPORTANT: The sub-graph has an additional node linked to sources (resp.
// destinations) for the forward (resp. backward) direction. This additional
// node is indexed with the last index. All added arcs are given to have an
// arc index in the original graph of -1.
const int sub_arcs_count = num_arcs + full_sources[dir].size();
sub_reverse_graph_[dir] = GraphType(sub_nodes.size() + 1, sub_arcs_count);
sub_arc_resources_[dir].resize(num_resources_);
for (int r = 0; r < num_resources_; ++r) {
sub_arc_resources_[dir][r].reserve(sub_arcs_count);
}
sub_full_arc_indices_[dir].reserve(sub_arcs_count);
const GraphType& reverse_full_graph = *(full_graph[Reverse(dir)]);
for (ArcIndex arc_index = 0; arc_index < num_arcs; ++arc_index) {
const NodeIndex from =
sub_node_indices_[dir][reverse_full_graph.Tail(arc_index)];
const NodeIndex to =
sub_node_indices_[dir][reverse_full_graph.Head(arc_index)];
if (from == -1 || to == -1) {
continue;
}
sub_reverse_graph_[dir].AddArc(from, to);
ArcIndex sub_full_arc_index;
if (dir == FORWARD && !full_inverse_arc_indices.empty()) {
sub_full_arc_index = full_inverse_arc_indices[arc_index];
} else {
sub_full_arc_index = arc_index;
}
for (int r = 0; r < num_resources_; ++r) {
sub_arc_resources_[dir][r].push_back(
(*arc_resources_)[r][sub_full_arc_index]);
}
sub_full_arc_indices_[dir].push_back(sub_full_arc_index);
}
for (const NodeIndex source : full_sources[dir]) {
const NodeIndex sub_source = sub_node_indices_[dir][source];
if (sub_source == -1) {
continue;
}
sub_reverse_graph_[dir].AddArc(sub_source, sub_nodes.size());
for (int r = 0; r < num_resources_; ++r) {
sub_arc_resources_[dir][r].push_back(0.0);
}
sub_full_arc_indices_[dir].push_back(-1);
}
std::vector<ArcIndex> sub_permutation;
sub_reverse_graph_[dir].Build(&sub_permutation);
for (int r = 0; r < num_resources_; ++r) {
util::Permute(sub_permutation, &sub_arc_resources_[dir][r]);
}
util::Permute(sub_permutation, &sub_full_arc_indices_[dir]);
// We add labels corresponding to the sources in the sub-graph as these will
// be artificially created labels due to the
max_num_created_labels_[dir] += full_sources[dir].size();
}
// Memory allocation is done here and only once in order to avoid
// reallocation at each call of `RunConstrainedShortestPathOnDag()` for
// better performance.
for (const Direction dir : {FORWARD, BACKWARD}) {
resources_from_sources_[dir].resize(num_resources_);
node_first_label_[dir].resize(sub_reverse_graph_[dir].size());
node_num_labels_[dir].resize(sub_reverse_graph_[dir].size());
}
}
template <class GraphType>
@@ -280,56 +520,308 @@ template <class GraphType>
#endif
PathWithLength ConstrainedShortestPathsOnDagWrapper<
GraphType>::RunConstrainedShortestPathOnDag() {
// Caching the vector addresses allow to not fetch it on each access.
const absl::Span<const double> arc_lengths = *arc_lengths_;
const absl::Span<const double> max_resources = *max_resources_;
// Clear all labels from previous run.
labels_from_sources_.clear();
double best_length = std::numeric_limits<double>::infinity();
int best_label_index = -1;
for (const NodeIndex to : *topological_order_) {
node_first_label_[to] = labels_from_sources_.size();
int& num_labels = node_num_labels_[to];
num_labels = 0;
if (node_is_after_last_destination_[to]) {
break;
// Assign lengths on sub-relevant graphs.
std::vector<double> sub_arc_lengths[2];
for (const Direction dir : {FORWARD, BACKWARD}) {
sub_arc_lengths[dir].reserve(sub_reverse_graph_[dir].num_arcs());
for (ArcIndex sub_arc_index = 0;
sub_arc_index < sub_reverse_graph_[dir].num_arcs(); ++sub_arc_index) {
const ArcIndex arc_index = sub_full_arc_indices_[dir][sub_arc_index];
if (arc_index == -1) {
sub_arc_lengths[dir].push_back(0.0);
continue;
}
sub_arc_lengths[dir].push_back((*arc_lengths_)[arc_index]);
}
if (node_is_source_[to]) {
labels_from_sources_.push_back(
{.length = 0,
.resources = std::vector<double>(max_resources.size()),
.incoming_arc = -1});
++num_labels;
if (labels_from_sources_.size() >= max_num_created_labels_) {
const std::vector<ArcIndex> arc_path =
BestArcPathEndingWith(best_label_index);
return {.length = best_length,
.arc_path = arc_path,
.node_path = NodePathImpliedBy(arc_path)};
}
// TODO(b/329036127): Each search could be done in parallel as no data is
// shared except `best_label_pair_`.
LabelPair best_label_pair = {
.length = std::numeric_limits<double>::infinity(),
.label_index = {-1, -1}};
for (const Direction dir : {FORWARD, BACKWARD}) {
RunHalfConstrainedShortestPathOnDag(
/*reverse_graph=*/sub_reverse_graph_[dir],
/*arc_lengths=*/sub_arc_lengths[dir],
/*arc_resources=*/sub_arc_resources_[dir],
/*min_arc_resources=*/sub_min_arc_resources_[dir],
/*max_resources=*/*max_resources_,
/*max_num_created_labels=*/max_num_created_labels_[dir],
/*lengths_from_sources=*/lengths_from_sources_[dir],
/*resources_from_sources=*/resources_from_sources_[dir],
/*incoming_arc_indices_from_sources=*/
incoming_arc_indices_from_sources_[dir],
/*first_label=*/node_first_label_[dir],
/*num_labels=*/node_num_labels_[dir]);
}
// Check destinations within relevant half sub-graphs.
for (const Direction dir : {FORWARD, BACKWARD}) {
absl::Span<const NodeIndex> destinations =
dir == FORWARD ? destinations_ : sources_;
for (const NodeIndex dst : destinations) {
const NodeIndex sub_dst = sub_node_indices_[dir][dst];
if (sub_dst == -1) {
continue;
}
const int num_labels_dst = node_num_labels_[dir][sub_dst];
if (num_labels_dst == 0) {
continue;
}
const int first_label_dst = node_first_label_[dir][sub_dst];
for (int label_index = first_label_dst;
label_index < first_label_dst + num_labels_dst; ++label_index) {
const double length_dst = lengths_from_sources_[dir][label_index];
if (length_dst < best_label_pair.length) {
best_label_pair.length = length_dst;
best_label_pair.label_index[dir] = label_index;
}
}
}
for (const ArcIndex reverse_arc_index : reverse_graph_.OutgoingArcs(to)) {
const NodeIndex from = reverse_graph_.Head(reverse_arc_index);
const ArcIndex arc_index =
reverse_inverse_arc_permutation_.empty()
? reverse_arc_index
: reverse_inverse_arc_permutation_[reverse_arc_index];
const double arc_length = arc_lengths[arc_index];
}
const ArcIndex merging_arc_index = MergeHalfRuns(
/*graph=*/*graph_, /*arc_lengths=*/*arc_lengths_,
/*arc_resources=*/*arc_resources_,
/*max_resources=*/*max_resources_,
/*sub_node_indices=*/sub_node_indices_,
/*lengths_from_sources=*/lengths_from_sources_,
/*resources_from_sources=*/resources_from_sources_,
/*first_label=*/node_first_label_,
/*num_labels=*/node_num_labels_, /*best_label_pair=*/best_label_pair);
std::vector<ArcIndex> arc_path;
for (const Direction dir : {FORWARD, BACKWARD}) {
for (const ArcIndex sub_arc_index : ArcPathTo(
/*best_label_index=*/best_label_pair.label_index[dir],
/*reverse_graph=*/sub_reverse_graph_[dir],
/*arc_lengths=*/sub_arc_lengths[dir],
/*lengths_from_sources=*/lengths_from_sources_[dir],
/*incoming_arc_indices_from_sources=*/
incoming_arc_indices_from_sources_[dir],
/*first_label=*/node_first_label_[dir],
/*num_labels=*/node_num_labels_[dir])) {
const ArcIndex arc_index = sub_full_arc_indices_[dir][sub_arc_index];
if (arc_index == -1) {
break;
}
arc_path.push_back(arc_index);
}
if (dir == FORWARD && merging_arc_index != -1) {
absl::c_reverse(arc_path);
arc_path.push_back(merging_arc_index);
}
}
// Clear all labels from the next run.
for (const Direction dir : {FORWARD, BACKWARD}) {
lengths_from_sources_[dir].clear();
for (int r = 0; r < num_resources_; ++r) {
resources_from_sources_[dir][r].clear();
}
incoming_arc_indices_from_sources_[dir].clear();
}
return {.length = best_label_pair.length,
.arc_path = arc_path,
.node_path = NodePathImpliedBy(arc_path, *graph_)};
}
template <class GraphType>
#if __cplusplus >= 202002L
requires DagGraphType<GraphType>
#endif
void ConstrainedShortestPathsOnDagWrapper<GraphType>::
RunHalfConstrainedShortestPathOnDag(
const GraphType& reverse_graph, absl::Span<const double> arc_lengths,
absl::Span<const std::vector<double>> arc_resources,
absl::Span<const std::vector<double>> min_arc_resources,
absl::Span<const double> max_resources,
const int max_num_created_labels,
std::vector<double>& lengths_from_sources,
std::vector<std::vector<double>>& resources_from_sources,
std::vector<ArcIndex>& incoming_arc_indices_from_sources,
std::vector<int>& first_label, std::vector<int>& num_labels) {
// Initialize source node.
const NodeIndex source_node = reverse_graph.num_nodes() - 1;
first_label[source_node] = 0;
num_labels[source_node] = 1;
lengths_from_sources.push_back(0);
for (int r = 0; r < num_resources_; ++r) {
resources_from_sources[r].push_back(0);
}
incoming_arc_indices_from_sources.push_back(-1);
std::vector<double> lengths_to;
std::vector<std::vector<double>> resources_to(num_resources_);
std::vector<ArcIndex> incoming_arc_indices_to;
std::vector<int> label_indices_to;
std::vector<double> resources(num_resources_);
for (NodeIndex to = 0; to < source_node; ++to) {
lengths_to.clear();
for (int r = 0; r < num_resources_; ++r) {
resources_to[r].clear();
}
incoming_arc_indices_to.clear();
for (const ArcIndex reverse_arc_index : reverse_graph.OutgoingArcs(to)) {
const NodeIndex from = reverse_graph.Head(reverse_arc_index);
const double arc_length = arc_lengths[reverse_arc_index];
DCHECK(arc_length != -std::numeric_limits<double>::infinity());
if (arc_length == std::numeric_limits<double>::infinity()) {
continue;
}
for (int i = node_first_label_[from];
i < node_first_label_[from] + node_num_labels_[from]; ++i) {
const Label& label_from = labels_from_sources_[i];
for (int label_index = first_label[from];
label_index < first_label[from] + num_labels[from]; ++label_index) {
bool path_is_feasible = true;
for (int r = 0; r < max_resources.size(); ++r) {
DCHECK_GE((*arc_resources_)[r][arc_index], 0.0);
// TODO(b/314756645): Include lower bound on the resources to
// increase pruning.
if (label_from.resources[r] + (*arc_resources_)[r][arc_index] >
for (int r = 0; r < num_resources_; ++r) {
DCHECK_GE(arc_resources[r][reverse_arc_index], 0.0);
resources[r] = resources_from_sources[r][label_index] +
arc_resources[r][reverse_arc_index];
if (resources[r] + min_arc_resources[r][to] > max_resources[r]) {
path_is_feasible = false;
break;
}
}
if (!path_is_feasible) {
continue;
}
lengths_to.push_back(lengths_from_sources[label_index] + arc_length);
for (int r = 0; r < num_resources_; ++r) {
resources_to[r].push_back(resources[r]);
}
incoming_arc_indices_to.push_back(reverse_arc_index);
}
}
// Sort labels lexicographically with lengths then resources.
label_indices_to.clear();
label_indices_to.reserve(lengths_to.size());
for (int i = 0; i < lengths_to.size(); ++i) {
label_indices_to.push_back(i);
}
absl::c_sort(label_indices_to, [&](const int i, const int j) {
if (lengths_to[i] < lengths_to[j]) return true;
if (lengths_to[i] > lengths_to[j]) return false;
for (int r = 0; r < num_resources_; ++r) {
if (resources_to[r][i] < resources_to[r][j]) return true;
if (resources_to[r][i] > resources_to[r][j]) return false;
}
return i < j;
});
first_label[to] = lengths_from_sources.size();
int& num_labels_to = num_labels[to];
// Reset the number of labels to zero otherwise it holds the previous run
// result.
num_labels_to = 0;
for (int i = 0; i < label_indices_to.size(); ++i) {
// Check if label "i" on node `to` is dominated by any other label.
const int label_i_index = label_indices_to[i];
bool label_i_is_dominated = false;
for (int j = 0; j < i - 1; ++j) {
const int label_j_index = label_indices_to[j];
if (lengths_to[label_i_index] <= lengths_to[label_j_index]) continue;
bool label_j_dominates_label_i = true;
for (int r = 0; r < num_resources_; ++r) {
if (resources_to[r][label_i_index] <=
resources_to[r][label_j_index]) {
label_j_dominates_label_i = false;
break;
}
}
if (label_j_dominates_label_i) {
label_i_is_dominated = true;
break;
}
}
if (label_i_is_dominated) continue;
lengths_from_sources.push_back(lengths_to[label_i_index]);
for (int r = 0; r < num_resources_; ++r) {
resources_from_sources[r].push_back(resources_to[r][label_i_index]);
}
incoming_arc_indices_from_sources.push_back(
incoming_arc_indices_to[label_i_index]);
++num_labels_to;
if (lengths_from_sources.size() >= max_num_created_labels) {
return;
}
}
}
}
template <class GraphType>
#if __cplusplus >= 202002L
requires DagGraphType<GraphType>
#endif
typename GraphType::ArcIndex
ConstrainedShortestPathsOnDagWrapper<GraphType>::MergeHalfRuns(
const GraphType& graph, absl::Span<const double> arc_lengths,
absl::Span<const std::vector<double>> arc_resources,
absl::Span<const double> max_resources,
const std::vector<NodeIndex> sub_node_indices[2],
const std::vector<double> lengths_from_sources[2],
const std::vector<std::vector<double>> resources_from_sources[2],
const std::vector<int> first_label[2], const std::vector<int> num_labels[2],
LabelPair& best_label_pair) {
const std::vector<NodeIndex>& forward_sub_node_indices =
sub_node_indices[FORWARD];
absl::Span<const double> forward_lengths = lengths_from_sources[FORWARD];
const std::vector<std::vector<double>>& forward_resources =
resources_from_sources[FORWARD];
absl::Span<const int> forward_first_label = first_label[FORWARD];
absl::Span<const int> forward_num_labels = num_labels[FORWARD];
const std::vector<NodeIndex>& backward_sub_node_indices =
sub_node_indices[BACKWARD];
absl::Span<const double> backward_lengths = lengths_from_sources[BACKWARD];
const std::vector<std::vector<double>>& backward_resources =
resources_from_sources[BACKWARD];
absl::Span<const int> backward_first_label = first_label[BACKWARD];
absl::Span<const int> backward_num_labels = num_labels[BACKWARD];
ArcIndex merging_arc_index = -1;
// TODO(b/328745703): These operations could be run in parallel with a mutex
// over `best_label_pair`.
for (ArcIndex arc_index = 0; arc_index < graph.num_arcs(); ++arc_index) {
const NodeIndex sub_from = forward_sub_node_indices[graph.Tail(arc_index)];
if (sub_from == -1) {
continue;
}
const NodeIndex sub_to = backward_sub_node_indices[graph.Head(arc_index)];
if (sub_to == -1) {
continue;
}
const int num_labels_from = forward_num_labels[sub_from];
if (num_labels_from == 0) {
continue;
}
const int num_labels_to = backward_num_labels[sub_to];
if (num_labels_to == 0) {
continue;
}
const double arc_length = arc_lengths[arc_index];
DCHECK(arc_length != -std::numeric_limits<double>::infinity());
if (arc_length == std::numeric_limits<double>::infinity()) {
continue;
}
const int first_label_from = forward_first_label[sub_from];
const int first_label_to = backward_first_label[sub_to];
for (int label_to_index = first_label_to;
label_to_index < first_label_to + num_labels_to; ++label_to_index) {
const double length_to = backward_lengths[label_to_index];
if (arc_length + length_to >= best_label_pair.length) {
continue;
}
for (int label_from_index = first_label_from;
label_from_index < first_label_from + num_labels_from;
++label_from_index) {
const double length_from = forward_lengths[label_from_index];
if (length_from + arc_length + length_to >= best_label_pair.length) {
continue;
}
bool path_is_feasible = true;
for (int r = 0; r < num_resources_; ++r) {
DCHECK_GE(arc_resources[r][arc_index], 0.0);
if (forward_resources[r][label_from_index] +
arc_resources[r][arc_index] +
backward_resources[r][label_to_index] >
max_resources[r]) {
path_is_feasible = false;
break;
@@ -338,67 +830,53 @@ PathWithLength ConstrainedShortestPathsOnDagWrapper<
if (!path_is_feasible) {
continue;
}
Label label_to = label_from;
label_to.length += arc_length;
for (int r = 0; r < max_resources.size(); ++r) {
label_to.resources[r] += (*arc_resources_)[r][arc_index];
}
label_to.incoming_arc = arc_index;
labels_from_sources_.push_back(label_to);
++num_labels;
if (node_is_destination_[to]) {
if (best_length > label_to.length) {
best_length = label_to.length;
best_label_index = labels_from_sources_.size() - 1;
}
}
if (labels_from_sources_.size() >= max_num_created_labels_) {
const std::vector<ArcIndex> arc_path =
BestArcPathEndingWith(best_label_index);
return {.length = best_length,
.arc_path = arc_path,
.node_path = NodePathImpliedBy(arc_path)};
}
best_label_pair.length = length_from + arc_length + length_to;
best_label_pair.label_index[FORWARD] = label_from_index;
best_label_pair.label_index[BACKWARD] = label_to_index;
merging_arc_index = arc_index;
}
}
}
const std::vector<ArcIndex> arc_path =
BestArcPathEndingWith(best_label_index);
return {.length = best_length,
.arc_path = arc_path,
.node_path = NodePathImpliedBy(arc_path)};
return merging_arc_index;
}
template <typename GraphType>
template <class GraphType>
#if __cplusplus >= 202002L
requires DagGraphType<GraphType>
#endif
std::vector<typename GraphType::ArcIndex>
ConstrainedShortestPathsOnDagWrapper<GraphType>::BestArcPathEndingWith(
const int label_index) const {
if (label_index < 0 || label_index >= labels_from_sources_.size()) {
ConstrainedShortestPathsOnDagWrapper<GraphType>::ArcPathTo(
const int best_label_index, const GraphType& reverse_graph,
absl::Span<const double> arc_lengths,
absl::Span<const double> lengths_from_sources,
absl::Span<const ArcIndex> incoming_arc_indices_from_sources,
absl::Span<const int> first_label, absl::Span<const int> num_labels) const {
if (best_label_index == -1) {
return {};
}
int current_label_index = label_index;
int current_label_index = best_label_index;
std::vector<ArcIndex> arc_path;
while (labels_from_sources_[current_label_index].incoming_arc != -1) {
const Label& label = labels_from_sources_[current_label_index];
const NodeIndex node = graph_->Tail(label.incoming_arc);
arc_path.push_back(label.incoming_arc);
for (int i = node_first_label_[node];
i < node_first_label_[node] + node_num_labels_[node]; ++i) {
// Since all labels of `labels_from_sources_` are feasible, we can pick
// any previous label which satisfies the length equality (irrespective of
// the resources it consumes).
if (labels_from_sources_[i].length +
(*arc_lengths_)[label.incoming_arc] ==
label.length) {
current_label_index = i;
for (int i = 0; i < reverse_graph.num_nodes(); ++i) {
const ArcIndex current_arc_index =
incoming_arc_indices_from_sources[current_label_index];
if (current_arc_index == -1) {
break;
}
arc_path.push_back(current_arc_index);
const NodeIndex sub_node = reverse_graph.Head(current_arc_index);
const double current_length = lengths_from_sources[current_label_index];
for (int label_index = first_label[sub_node];
label_index < first_label[sub_node] + num_labels[sub_node];
++label_index) {
if (std::abs(lengths_from_sources[label_index] +
arc_lengths[current_arc_index] - current_length) <=
kTolerance) {
current_label_index = label_index;
break;
}
}
}
absl::c_reverse(arc_path);
CHECK_EQ(incoming_arc_indices_from_sources[current_label_index], -1);
return arc_path;
}
@@ -408,16 +886,16 @@ template <typename GraphType>
#endif
std::vector<typename GraphType::NodeIndex>
ConstrainedShortestPathsOnDagWrapper<GraphType>::NodePathImpliedBy(
const std::vector<ArcIndex>& arc_path) const {
absl::Span<const ArcIndex> arc_path, const GraphType& graph) const {
if (arc_path.empty()) {
return {};
}
std::vector<NodeIndex> node_path;
node_path.reserve(arc_path.size() + 1);
for (const ArcIndex arc_index : arc_path) {
node_path.push_back(graph_->Tail(arc_index));
node_path.push_back(graph.Tail(arc_index));
}
node_path.push_back(graph_->Head(arc_path.back()));
node_path.push_back(graph.Head(arc_path.back()));
return node_path;
}

194
ortools/graph/dag_constrained_shortest_path_test.cc
View File

@@ -22,6 +22,7 @@
#include "absl/algorithm/container.h"
#include "absl/log/check.h"
#include "absl/random/random.h"
#include "absl/strings/str_cat.h"
#include "absl/types/span.h"
#include "benchmark/benchmark.h"
#include "gtest/gtest.h"
@@ -30,6 +31,7 @@
#include "ortools/graph/dag_shortest_path.h"
#include "ortools/graph/graph.h"
#include "ortools/graph/io.h"
#include "ortools/math_opt/cpp/math_opt.h"
namespace operations_research {
namespace {
@@ -154,6 +156,21 @@ TEST(ConstrainedShortestPathOnDagTest, GraphWithInefficientEdge) {
/*node_path=*/ElementsAre(source, a, destination)));
}
TEST(ConstrainedShortestPathOnDagTest, NoResources) {
const int source = 0;
const int destination = 1;
const int a = 2;
const int num_nodes = 3;
const std::vector<ArcWithLengthAndResources> arcs_with_length_and_resources =
{{source, a, 5.0, {}}, {a, destination, 3.0, {}}};
EXPECT_DEATH(
ConstrainedShortestPathsOnDag(num_nodes, arcs_with_length_and_resources,
source, destination,
/*max_resources=*/{}),
"ortools/graph/dag_shortest_path.h");
}
TEST(ConstrainedShortestPathOnDagTest, Cycle) {
const int source = 0;
const int destination = 1;
@@ -278,19 +295,21 @@ TEST(ConstrainedShortestPathsOnDagWrapperTest,
const int num_nodes = 3;
util::ListGraph<> graph(num_nodes, /*arc_capacity=*/2);
std::vector<double> arc_lengths;
std::vector<std::vector<double>> arc_resources;
std::vector<std::vector<double>> arc_resources(1);
graph.AddArc(source_1, source_2);
arc_lengths.push_back(-7.0);
arc_resources[0].push_back(2.0);
graph.AddArc(source_2, destination);
arc_lengths.push_back(3.0);
arc_resources[0].push_back(3.0);
const std::vector<int> topological_order = {source_1, source_2, destination};
const std::vector<int> sources = {source_1, source_2};
const std::vector<int> destinations = {destination};
const std::vector<double> max_resources = {};
const std::vector<double> max_resources = {6.0};
ConstrainedShortestPathsOnDagWrapper<util::ListGraph<>>
constrained_shortest_path_on_dag(&graph, &arc_lengths, &arc_resources,
&topological_order, &sources,
&destinations, &max_resources);
topological_order, sources, destinations,
&max_resources);
EXPECT_THAT(
constrained_shortest_path_on_dag.RunConstrainedShortestPathOnDag(),
@@ -324,8 +343,8 @@ TEST(ConstrainedShortestPathsOnDagWrapperTest, MultipleDestinations) {
const std::vector<double> max_resources = {6.0};
ConstrainedShortestPathsOnDagWrapper<util::ListGraph<>>
constrained_shortest_path_on_dag(&graph, &arc_lengths, &arc_resources,
&topological_order, &sources,
&destinations, &max_resources);
topological_order, sources, destinations,
&max_resources);
EXPECT_THAT(
constrained_shortest_path_on_dag.RunConstrainedShortestPathOnDag(),
@@ -364,8 +383,8 @@ TEST(ConstrainedShortestPathsOnDagWrapperTest, UpdateArcsLength) {
const std::vector<double> max_resources = {6.0, 12.0};
ConstrainedShortestPathsOnDagWrapper<util::ListGraph<>>
constrained_shortest_path_on_dag(&graph, &arc_lengths, &arc_resources,
&topological_order, &sources,
&destinations, &max_resources);
topological_order, sources, destinations,
&max_resources);
EXPECT_THAT(
constrained_shortest_path_on_dag.RunConstrainedShortestPathOnDag(),
@@ -397,30 +416,22 @@ TEST(ConstrainedShortestPathsOnDagWrapperTest, LimitMaximumNumberOfLabels) {
arc_resources[0].push_back(1.0);
const std::vector<int> topological_order = {source, a, destination};
const std::vector<int> sources = {source};
const std::vector<int> destinations = {a, destination};
const std::vector<int> destinations = {destination};
const std::vector<double> max_resources = {6.0};
ConstrainedShortestPathsOnDagWrapper<util::ListGraph<>>
constrained_shortest_path_on_dag_with_one_label(
&graph, &arc_lengths, &arc_resources, &topological_order, &sources,
&destinations, &max_resources, /*max_num_created_labels=*/1);
ConstrainedShortestPathsOnDagWrapper<util::ListGraph<>>
constrained_shortest_path_on_dag_with_two_labels(
&graph, &arc_lengths, &arc_resources, &topological_order, &sources,
&destinations, &max_resources, /*max_num_created_labels=*/2);
&graph, &arc_lengths, &arc_resources, topological_order, sources,
destinations, &max_resources, /*max_num_created_labels=*/1);
ConstrainedShortestPathsOnDagWrapper<util::ListGraph<>>
constrained_shortest_path_on_dag_with_three_labels(
&graph, &arc_lengths, &arc_resources, &topological_order, &sources,
&destinations, &max_resources, /*max_num_created_labels=*/3);
&graph, &arc_lengths, &arc_resources, topological_order, sources,
destinations, &max_resources, /*max_num_created_labels=*/3);
EXPECT_THAT(constrained_shortest_path_on_dag_with_one_label
.RunConstrainedShortestPathOnDag(),
FieldsAre(/*length=*/kInf,
/*arc_path=*/IsEmpty(),
/*node_path=*/IsEmpty()));
EXPECT_THAT(constrained_shortest_path_on_dag_with_two_labels
.RunConstrainedShortestPathOnDag(),
FieldsAre(/*length=*/5.0, /*arc_path=*/ElementsAre(0),
/*node_path=*/ElementsAre(source, a)));
EXPECT_THAT(constrained_shortest_path_on_dag_with_three_labels
.RunConstrainedShortestPathOnDag(),
FieldsAre(/*length=*/3.0, /*arc_path=*/ElementsAre(0, 1),
@@ -478,11 +489,67 @@ std::vector<double> GenerateRandomIntegerValues(
return arc_values;
}
// TODO(b/316203403): Remplace bounds with correct value computed using MIP.
double SolveConstrainedShortestPathUsingIntegerProgramming(
const util::StaticGraph<>& graph, absl::Span<const double> arc_lengths,
absl::Span<const std::vector<double>> arc_resources,
absl::Span<const double> max_resources,
absl::Span<const util::StaticGraph<>::NodeIndex> sources,
absl::Span<const util::StaticGraph<>::NodeIndex> destinations) {
using NodeIndex = util::StaticGraph<>::NodeIndex;
using ArcIndex = util::StaticGraph<>::ArcIndex;
math_opt::Model model;
std::vector<math_opt::Variable> arc_variables;
std::vector<math_opt::LinearExpression> flow_conservation(graph.num_nodes(),
0.0);
for (ArcIndex arc_index = 0; arc_index < graph.num_arcs(); ++arc_index) {
arc_variables.push_back(model.AddBinaryVariable(absl::StrCat(
arc_index, "_", graph.Tail(arc_index), "->", graph.Head(arc_index))));
model.set_objective_coefficient(arc_variables[arc_index],
arc_lengths[arc_index]);
flow_conservation[graph.Head(arc_index)] -= arc_variables[arc_index];
flow_conservation[graph.Tail(arc_index)] += arc_variables[arc_index];
}
math_opt::LinearExpression all_sources;
math_opt::LinearExpression all_destinations;
for (NodeIndex node_index = 0; node_index < graph.num_nodes(); ++node_index) {
math_opt::LinearExpression net_flow = 0;
if (absl::c_linear_search(sources, node_index)) {
const math_opt::Variable s = model.AddBinaryVariable();
all_sources += s;
net_flow += s;
}
if (absl::c_linear_search(destinations, node_index)) {
const math_opt::Variable t = model.AddBinaryVariable();
all_destinations += t;
net_flow -= t;
}
model.AddLinearConstraint(flow_conservation[node_index] == net_flow);
}
model.AddLinearConstraint(all_sources == 1);
model.AddLinearConstraint(all_destinations == 1);
for (int r = 0; r < max_resources.size(); ++r) {
math_opt::LinearExpression variable_resources;
for (ArcIndex arc_index = 0; arc_index < graph.num_arcs(); ++arc_index) {
variable_resources +=
arc_resources[r][arc_index] * arc_variables[arc_index];
}
model.AddLinearConstraint(variable_resources <= max_resources[r]);
}
const absl::StatusOr<math_opt::SolveResult> result =
math_opt::Solve(model, math_opt::SolverType::kCpSat, {});
CHECK_OK(result.status())
<< util::GraphToString(graph, util::PRINT_GRAPH_ARCS);
CHECK_OK(result->termination.EnsureIsOptimal())
<< util::GraphToString(graph, util::PRINT_GRAPH_ARCS);
return result->objective_value();
}
TEST(ConstrainedShortestPathsOnDagWrapperTest,
RandomizedStressTestSingleResource) {
absl::BitGen bit_gen;
const int kNumTests = 10000;
const int kNumTests = 50;
for (int test = 0; test < kNumTests; ++test) {
const int num_nodes = absl::Uniform(bit_gen, 2, 12);
const int num_arcs = absl::Uniform(
@@ -495,55 +562,23 @@ TEST(ConstrainedShortestPathsOnDagWrapperTest,
/*max_value=*/10.0,
/*start_to_end_value=*/1.0);
const std::vector<int> sources = {0};
std::vector<int> destinations(num_nodes);
absl::c_iota(destinations, 0);
const std::vector<int> destinations = {num_nodes - 1};
const std::vector<double> max_resources = {15.0};
ConstrainedShortestPathsOnDagWrapper<util::StaticGraph<>>
constrained_shortest_path_on_dag(&graph, &arc_lengths, &arc_resources,
&topological_order, &sources,
&destinations, &max_resources);
// Run DAG shortest path on `arc_length` as a LB.
ShortestPathsOnDagWrapper<util::StaticGraph<>> shortest_path_length(
&graph, &arc_lengths, topological_order);
// Run DAG shortest path on `arc_resources` as a UB.
ShortestPathsOnDagWrapper<util::StaticGraph<>> shortest_path_resources(
&graph, &(arc_resources[0]), topological_order);
shortest_path_resources.RunShortestPathOnDag(sources);
// Precompute the expected reached nodes: any node whose resource is <=
// max_resources[0].
std::vector<int> expected_reached_nodes;
for (int node = 0; node < num_nodes; ++node) {
if (shortest_path_resources.LengthTo(node) <= max_resources[0]) {
expected_reached_nodes.push_back(node);
}
}
const int kNumSamples = 20;
topological_order, sources,
destinations, &max_resources);
const int kNumSamples = 5;
for (int _ = 0; _ < kNumSamples; ++_) {
// Draw random lengths and recompute *un*constrained shortest paths.
arc_lengths = GenerateRandomIntegerValues(graph);
shortest_path_length.RunShortestPathOnDag(sources);
const PathWithLength path_with_length =
constrained_shortest_path_on_dag.RunConstrainedShortestPathOnDag();
EXPECT_FALSE(path_with_length.node_path.empty());
const int best_destination = path_with_length.node_path.back();
const std::vector<int> ub_arc_path =
shortest_path_resources.ArcPathTo(best_destination);
double ub_length = 0.0;
for (const int arc : ub_arc_path) {
ub_length += arc_lengths[arc];
}
EXPECT_GE(path_with_length.length,
shortest_path_length.LengthTo(best_destination))
<< best_destination;
EXPECT_LE(path_with_length.length, ub_length) << best_destination;
EXPECT_NEAR(path_with_length.length,
SolveConstrainedShortestPathUsingIntegerProgramming(
graph, arc_lengths, arc_resources, max_resources, sources,
destinations),
1e-5);
ASSERT_FALSE(HasFailure())
<< DUMP_VARS(num_nodes, num_arcs, arc_lengths) << "\n With graph :\n "
@@ -578,25 +613,33 @@ void BM_RandomDag(benchmark::State& state) {
const std::vector<double> max_resources = {num_nodes * 0.2};
ConstrainedShortestPathsOnDagWrapper<util::StaticGraph<>>
constrained_shortest_path_on_dag(&graph, &arc_lengths, &arc_resources,
&topological_order, &sources,
&destinations, &max_resources);
topological_order, sources, destinations,
&max_resources);
int total_label_count = 0;
for (auto _ : state) {
// Pick a arc lengths scenario at random.
arc_lengths =
arc_lengths_scenarios[absl::Uniform(bit_gen, 0, num_scenarios)];
const PathWithLength path_with_length =
constrained_shortest_path_on_dag.RunConstrainedShortestPathOnDag();
total_label_count += constrained_shortest_path_on_dag.label_count();
CHECK_GE(path_with_length.length, 0.0);
CHECK_LE(path_with_length.length, 10000.0);
}
state.SetItemsProcessed(state.iterations() * (num_nodes + num_arcs));
state.SetItemsProcessed(std::max(1, total_label_count));
}
BENCHMARK(BM_RandomDag)
->ArgPair(1 << 6, 4)
->ArgPair(1 << 6, 16)
->ArgPair(1 << 9, 4);
->ArgPair(1 << 9, 4)
->ArgPair(1 << 9, 16)
->ArgPair(1 << 12, 4)
->ArgPair(1 << 12, 16)
->ArgPair(1 << 15, 4)
->ArgPair(1 << 15, 16)
->ArgPair(1 << 18, 4);
// -----------------------------------------------------------------------------
// Debug tests.
@@ -654,6 +697,17 @@ TEST(ConstrainedShortestPathOnDagTest, NegativeMaxResource) {
"negative");
}
TEST(ConstrainedShortestPathOnDagTest, SourceIsDestination) {
const int source = 0;
const int num_nodes = 1;
EXPECT_DEATH(
ConstrainedShortestPathsOnDag(
num_nodes, /*arcs_with_length_and_resources=*/{}, source, source,
/*max_resources=*/{0.0}),
"source and destination");
}
TEST(ConstrainedShortestPathsOnDagWrapperTest, ValidateTopologicalOrder) {
const int source = 0;
const int destination = 1;
@@ -670,8 +724,8 @@ TEST(ConstrainedShortestPathsOnDagWrapperTest, ValidateTopologicalOrder) {
const std::vector<double> max_resources = {0.0};
EXPECT_DEATH(ConstrainedShortestPathsOnDagWrapper<util::ListGraph<>>(
&graph, &arc_lengths, &arc_resources, &topological_order,
&sources, &destinations, &max_resources),
&graph, &arc_lengths, &arc_resources, topological_order,
sources, destinations, &max_resources),
"Invalid topological order");
}
#endif // NDEBUG