From 92be18748bd93102e3c4e707345bdf657f08c1de Mon Sep 17 00:00:00 2001 From: Laurent Perron Date: Mon, 20 Jan 2025 14:58:09 +0100 Subject: [PATCH] more graph cleanup --- examples/cpp/max_flow.cc | 2 +- ortools/graph/generic_max_flow.h | 270 ++++++++++++++----------- ortools/graph/generic_max_flow_test.cc | 80 +++++++- ortools/graph/shortest_paths.h | 8 - 4 files changed, 232 insertions(+), 128 deletions(-) diff --git a/examples/cpp/max_flow.cc b/examples/cpp/max_flow.cc index c4171d4ff2..f4b4f18099 100644 --- a/examples/cpp/max_flow.cc +++ b/examples/cpp/max_flow.cc @@ -26,7 +26,7 @@ using Graph = ::util::ReverseArcListGraph<>; using NodeIndex = Graph::NodeIndex; using ArcIndex = Graph::ArcIndex; using MaxFlowT = GenericMaxFlow; -using FlowQuantity = MaxFlowT::FlowQuantityT; +using FlowQuantity = MaxFlowT::FlowSumType; void SolveMaxFlow() { const int num_nodes = 5; diff --git a/ortools/graph/generic_max_flow.h b/ortools/graph/generic_max_flow.h index bd068ed29d..c03a964b9c 100644 --- a/ortools/graph/generic_max_flow.h +++ b/ortools/graph/generic_max_flow.h @@ -230,12 +230,31 @@ class MaxFlowStatusClass { // reverse can have a positive initial capacity. This can lead to twice fewer // arcs and a faster algo if the "user graph" as a lot of reverse arcs // already. -template +// +// Flow types and overflow: To optimize memory usage and speed we have two +// integer types representing flow quantities. +// - ArcFlowType (can be unsigned) is used to represent the flow per arc. The +// type should be big enough to hold the initial capacity of an arc + its +// reverse. So if all your initial capacity <= +// std::numeric_limit::max() / 2, you are safe. +// - FlowSumType (must be signed) is used to hold the sum of a flow going +// through a node or through the full network. Ideally, the sum of all the +// capacity of the arcs leaving the source should fit in a FlowSumType. If +// this is the case, you will get the better performance, and always an +// OPTIMAL result. Otherwise, the code might be slower, and you will get +// OPTIMAL only if the max_flow that can be sent through the network fit a +// FlowSumType, you will get the INT_OVERFLOW status if more flow than that +// can be sent. +template class GenericMaxFlow : public MaxFlowStatusClass { public: typedef typename Graph::NodeIndex NodeIndex; typedef typename Graph::ArcIndex ArcIndex; - typedef FlowType FlowQuantityT; + typedef ArcFlowT ArcFlowType; + typedef FlowSumT FlowSumType; + static_assert(sizeof(FlowSumType) >= sizeof(ArcFlowType)); + static_assert(std::is_signed_v); // The height of a node never excess 2 times the number of node, so we // use the same type as a Node index. @@ -253,8 +272,6 @@ class GenericMaxFlow : public MaxFlowStatusClass { GenericMaxFlow& operator=(const GenericMaxFlow&) = delete; #endif - virtual ~GenericMaxFlow() {} - // Returns the graph associated to the current object. const Graph* graph() const { return graph_; } @@ -273,29 +290,36 @@ class GenericMaxFlow : public MaxFlowStatusClass { // // Note that this will be ignored for self-arc, so do not be surprised if you // get zero when reading the capacity of a self-arc back. - void SetArcCapacity(ArcIndex arc, FlowType new_capacity); + void SetArcCapacity(ArcIndex arc, ArcFlowType new_capacity); // Returns true if a maximum flow was solved. bool Solve(); // Returns the total flow found by the algorithm. - FlowType GetOptimalFlow() const { return node_excess_[sink_]; } + FlowSumType GetOptimalFlow() const { return node_excess_[sink_]; } // Returns the current flow on the given arc. // // Note that on a couple (arc, opposite_arc) the flow only goes in one // direction (where it is positive) and the other direction will have the // negation of that flow. - FlowType Flow(ArcIndex arc) const { + // + // Tricky: This returns a FlowSumType in order to be able to distinguish + // negative flow from positive flow when ArcFlowType is unsigned. + // TODO(user): We could change the API to avoid that if needed. + FlowSumType Flow(ArcIndex arc) const { if constexpr (Graph::kHasNegativeReverseArcs) { - if (IsArcDirect(arc)) return residual_arc_capacity_[Opposite(arc)]; - return -residual_arc_capacity_[arc]; + if (IsArcDirect(arc)) { + return static_cast(residual_arc_capacity_[Opposite(arc)]); + } + return -static_cast(residual_arc_capacity_[arc]); } - return initial_capacity_[arc] - residual_arc_capacity_[arc]; + return static_cast(initial_capacity_[arc]) - + static_cast(residual_arc_capacity_[arc]); } // Returns the initial capacity of an arc. - FlowType Capacity(ArcIndex arc) const { + ArcFlowType Capacity(ArcIndex arc) const { if constexpr (Graph::kHasNegativeReverseArcs) { if (!IsArcDirect(arc)) return 0; return residual_arc_capacity_[arc] + @@ -401,19 +425,22 @@ class GenericMaxFlow : public MaxFlowStatusClass { // Tries to saturate all the outgoing arcs from the source that can reach the // sink. Most of the time, we can do that in one go, except when more flow - // than kMaxFlowQuantity can be pushed out of the source in which case we have + // than kMaxFlowSum can be pushed out of the source in which case we have // to be careful. Returns true if some flow was pushed. bool SaturateOutgoingArcsFromSource(); // Pushes flow on arc, i.e. consumes flow on residual_arc_capacity_[arc], // and consumes -flow on residual_arc_capacity_[Opposite(arc)]. Updates // node_excess_ at the tail and head of arc accordingly. - void PushFlow(FlowType flow, NodeIndex tail, ArcIndex arc); + void PushFlow(ArcFlowType flow, NodeIndex tail, ArcIndex arc); - // Same as PushFlow() but with a cached node_excess_.data(), this is faster - // in hot loop as we remove bound checking and the pointer is constant. - void FastPushFlow(FlowType flow, NodeIndex tail, ArcIndex arc, - FlowType* node_excess); + // Similar to PushFlow(), but this will always push + // std::min(node_excess_[tail], residual_arc_capacity_[arc]). + // + // We also use a cached node_excess_.data() as this is used in hot loops and + // allow to remove bound checking and fetching the pointer which is constant. + void PushAsMuchFlowAsPossible(NodeIndex tail, ArcIndex arc, + FlowSumType* node_excess); // Relabels a node, i.e. increases its height by the minimum necessary amount. // This version of Relabel is relaxed in a way such that if an admissible arc @@ -435,14 +462,16 @@ class GenericMaxFlow : public MaxFlowStatusClass { void ComputeReachableNodes(NodeIndex start, std::vector* result); // Maximum manageable flow. - static constexpr FlowType kMaxFlowQuantity = - std::numeric_limits::max(); + static constexpr FlowSumType kMaxFlowSum = + std::numeric_limits::max(); // A pointer to the graph passed as argument. const Graph* graph_; // An array representing the excess for each node in graph_. - std::unique_ptr node_excess_; + // With the current algorithm this is always positive except at the source. + // TODO(user): If needed we could tweak that to allow unsigned FlowSumType. + std::unique_ptr node_excess_; // An array representing the height function for each node in graph_. For a // given node, this is a lower bound on the shortest path length from this @@ -473,12 +502,12 @@ class GenericMaxFlow : public MaxFlowStatusClass { // Using these facts enables one to only maintain residual_arc_capacity_, // instead of both capacity and flow, for each direct and indirect arc. This // reduces the amount of memory for this information by a factor 2. - ZVector residual_arc_capacity_; + ZVector residual_arc_capacity_; // The initial capacity as set by SetArcCapacity(), unused if // `Graph::kNegativeReverseArcs`, as we can always recover the initial // capacity from the residual capacities of an arc and its reverse. - std::unique_ptr initial_capacity_; + std::unique_ptr initial_capacity_; // An array representing the first admissible arc for each node in graph_. std::unique_ptr first_admissible_arc_; @@ -560,10 +589,10 @@ Element PriorityQueueWithRestrictedPush::PopBack( return element; } -template -GenericMaxFlow::GenericMaxFlow(const Graph* graph, - NodeIndex source, - NodeIndex sink) +template +GenericMaxFlow::GenericMaxFlow(const Graph* graph, + NodeIndex source, + NodeIndex sink) : graph_(graph), residual_arc_capacity_(), source_(source), @@ -578,7 +607,7 @@ GenericMaxFlow::GenericMaxFlow(const Graph* graph, // // TODO(user): Unfortunately std/absl::make_unique_for_overwrite is not // yet available in open-source. - node_excess_ = std::make_unique(max_num_nodes); + node_excess_ = std::make_unique(max_num_nodes); node_potential_ = std::make_unique(max_num_nodes); first_admissible_arc_ = std::make_unique(max_num_nodes); bfs_queue_.reserve(max_num_nodes); @@ -589,16 +618,16 @@ GenericMaxFlow::GenericMaxFlow(const Graph* graph, residual_arc_capacity_.Reserve(-max_num_arcs, max_num_arcs - 1); } else { // We will need to store the initial capacity in this case. - initial_capacity_ = std::make_unique(max_num_arcs); + initial_capacity_ = std::make_unique(max_num_arcs); residual_arc_capacity_.Reserve(0, max_num_arcs - 1); } residual_arc_capacity_.SetAll(0); } } -template -void GenericMaxFlow::SetArcCapacity(ArcIndex arc, - FlowType new_capacity) { +template +void GenericMaxFlow::SetArcCapacity( + ArcIndex arc, ArcFlowType new_capacity) { SCOPED_TIME_STAT(&stats_); DCHECK_LE(0, new_capacity); DCHECK(IsArcDirect(arc)); @@ -616,23 +645,25 @@ void GenericMaxFlow::SetArcCapacity(ArcIndex arc, residual_arc_capacity_[Opposite(arc)] = 0; } else { initial_capacity_[arc] = new_capacity; + DCHECK_LE(initial_capacity_[arc], std::numeric_limits::max() - + initial_capacity_[Opposite(arc)]); } } -template -void GenericMaxFlow::GetSourceSideMinCut( +template +void GenericMaxFlow::GetSourceSideMinCut( std::vector* result) { ComputeReachableNodes(source_, result); } -template -void GenericMaxFlow::GetSinkSideMinCut( +template +void GenericMaxFlow::GetSinkSideMinCut( std::vector* result) { ComputeReachableNodes(sink_, result); } -template -bool GenericMaxFlow::CheckResult() const { +template +bool GenericMaxFlow::CheckResult() const { SCOPED_TIME_STAT(&stats_); if (node_excess_[source_] != -node_excess_[sink_]) { LOG(DFATAL) << "-node_excess_[source_] = " << -node_excess_[source_] @@ -650,8 +681,8 @@ bool GenericMaxFlow::CheckResult() const { } for (ArcIndex arc = 0; arc < graph_->num_arcs(); ++arc) { const ArcIndex opposite = Opposite(arc); - const FlowType direct_capacity = residual_arc_capacity_[arc]; - const FlowType opposite_capacity = residual_arc_capacity_[opposite]; + const ArcFlowType direct_capacity = residual_arc_capacity_[arc]; + const ArcFlowType opposite_capacity = residual_arc_capacity_[opposite]; if (direct_capacity < 0) { LOG(DFATAL) << "residual_arc_capacity_[" << arc << "] = " << direct_capacity << " < 0"; @@ -670,7 +701,7 @@ bool GenericMaxFlow::CheckResult() const { } } - if (GetOptimalFlow() < kMaxFlowQuantity && AugmentingPathExists()) { + if (GetOptimalFlow() < kMaxFlowSum && AugmentingPathExists()) { LOG(DFATAL) << "The algorithm terminated, but the flow is not maximal!"; return false; } @@ -678,8 +709,8 @@ bool GenericMaxFlow::CheckResult() const { return true; } -template -bool GenericMaxFlow::AugmentingPathExists() const { +template +bool GenericMaxFlow::AugmentingPathExists() const { SCOPED_TIME_STAT(&stats_); // We simply compute the reachability from the source in the residual graph. @@ -705,8 +736,8 @@ bool GenericMaxFlow::AugmentingPathExists() const { return is_reached[sink_]; } -template -bool GenericMaxFlow::CheckRelabelPrecondition( +template +bool GenericMaxFlow::CheckRelabelPrecondition( NodeIndex node) const { DCHECK(IsActive(node)); for (const ArcIndex arc : graph_->OutgoingOrOppositeIncomingArcs(node)) { @@ -716,8 +747,8 @@ bool GenericMaxFlow::CheckRelabelPrecondition( return true; } -template -std::string GenericMaxFlow::DebugString( +template +std::string GenericMaxFlow::DebugString( absl::string_view context, ArcIndex arc) const { const NodeIndex tail = Tail(arc); const NodeIndex head = Head(arc); @@ -732,8 +763,8 @@ std::string GenericMaxFlow::DebugString( node_potential_[head], node_excess_[tail], node_excess_[head]); } -template -bool GenericMaxFlow::Solve() { +template +bool GenericMaxFlow::Solve() { status_ = NOT_SOLVED; InitializePreflow(); @@ -753,16 +784,16 @@ bool GenericMaxFlow::Solve() { status_ = OPTIMAL; DCHECK(CheckResult()); - if (GetOptimalFlow() == kMaxFlowQuantity && AugmentingPathExists()) { - // In this case, we are sure that the flow is > kMaxFlowQuantity. + if (GetOptimalFlow() == kMaxFlowSum && AugmentingPathExists()) { + // In this case, we are sure that the flow is > kMaxFlowSum. status_ = INT_OVERFLOW; } IF_STATS_ENABLED(VLOG(1) << stats_.StatString()); return true; } -template -void GenericMaxFlow::InitializePreflow() { +template +void GenericMaxFlow::InitializePreflow() { SCOPED_TIME_STAT(&stats_); // TODO(user): Ebert graph has an issue with nodes with no arcs, so we @@ -810,8 +841,8 @@ void GenericMaxFlow::InitializePreflow() { // potentials because of the way we cancel flow on cycle. However, we only call // that at the end of the algorithm, or just before a GlobalUpdate() that will // restore the precondition on the node potentials. -template -void GenericMaxFlow::PushFlowExcessBackToSource() { +template +void GenericMaxFlow::PushFlowExcessBackToSource() { SCOPED_TIME_STAT(&stats_); const NodeIndex num_nodes = graph_->num_nodes(); @@ -880,7 +911,7 @@ void GenericMaxFlow::PushFlowExcessBackToSource() { visited[node] = true; index_branch.push_back(arc_stack.size() - 1); for (const ArcIndex arc : graph_->OutgoingArcs(node)) { - const FlowType flow = Flow(arc); + const FlowSumType flow = Flow(arc); const NodeIndex head = Head(arc); if (flow > 0 && !stored[head]) { if (!visited[head]) { @@ -898,11 +929,11 @@ void GenericMaxFlow::PushFlowExcessBackToSource() { // Compute the maximum flow that can be canceled on the cycle and the // min index such that arc_stack[index_branch[i]] will be saturated. - FlowType flow_on_cycle = flow; + ArcFlowType flow_on_cycle = static_cast(flow); int first_saturated_index = index_branch.size(); for (int i = index_branch.size() - 1; i >= cycle_begin; --i) { const ArcIndex arc_on_cycle = arc_stack[index_branch[i]]; - if (const FlowType arc_flow = Flow(arc_on_cycle); + if (const ArcFlowType arc_flow = Flow(arc_on_cycle); arc_flow <= flow_on_cycle) { flow_on_cycle = arc_flow; first_saturated_index = i; @@ -910,7 +941,7 @@ void GenericMaxFlow::PushFlowExcessBackToSource() { } // This is just here for a DCHECK() below. - const FlowType excess = node_excess_[head]; + const FlowSumType excess = node_excess_[head]; // Cancel the flow on the cycle, and set visited[node] = false for // the node that will be backtracked over. @@ -952,10 +983,11 @@ void GenericMaxFlow::PushFlowExcessBackToSource() { const NodeIndex node = reverse_topological_order[i]; if (node_excess_[node] == 0) continue; for (const ArcIndex arc : graph_->OutgoingOrOppositeIncomingArcs(node)) { - const FlowType flow = Flow(arc); + const FlowSumType flow = Flow(arc); if (flow < 0) { DCHECK_GT(residual_arc_capacity_[arc], 0); - const FlowType to_push = std::min(node_excess_[node], -flow); + const ArcFlowType to_push = + static_cast(std::min(node_excess_[node], -flow)); PushFlow(to_push, node, arc); if (node_excess_[node] == 0) break; } @@ -965,8 +997,8 @@ void GenericMaxFlow::PushFlowExcessBackToSource() { DCHECK_EQ(-node_excess_[source_], node_excess_[sink_]); } -template -void GenericMaxFlow::GlobalUpdate() { +template +void GenericMaxFlow::GlobalUpdate() { SCOPED_TIME_STAT(&stats_); bfs_queue_.clear(); int queue_index = 0; @@ -975,7 +1007,7 @@ void GenericMaxFlow::GlobalUpdate() { node_in_bfs_queue_[sink_] = true; // We cache this as this is a hot-loop and we don't want bound checking. - FlowType* const node_excess = node_excess_.get(); + FlowSumType* const node_excess = node_excess_.get(); // We do a BFS in the reverse residual graph, starting from the sink. // Because all the arcs from the source are saturated (except in @@ -1006,7 +1038,7 @@ void GenericMaxFlow::GlobalUpdate() { if (residual_arc_capacity_[opposite_arc] > 0) { // Note(user): We used to have a DCHECK_GE(candidate_distance, // node_potential_[head]); which is always true except in the case - // where we can push more than kMaxFlowQuantity out of the source. The + // where we can push more than kMaxFlowSum out of the source. The // problem comes from the fact that in this case, we call // PushFlowExcessBackToSource() in the middle of the algorithm. The // later call will break the properties of the node potential. Note @@ -1018,9 +1050,7 @@ void GenericMaxFlow::GlobalUpdate() { // Note(user): I haven't seen this anywhere in the literature. // TODO(user): Investigate more and maybe write a publication :) if (node_excess[head] > 0) { - const FlowType flow = - std::min(node_excess[head], residual_arc_capacity_[opposite_arc]); - FastPushFlow(flow, head, opposite_arc, node_excess); + PushAsMuchFlowAsPossible(head, opposite_arc, node_excess); // If the arc became saturated, it is no longer in the residual // graph, so we do not need to consider head at this time. @@ -1065,38 +1095,39 @@ void GenericMaxFlow::GlobalUpdate() { } } -template -bool GenericMaxFlow::SaturateOutgoingArcsFromSource() { +template +bool GenericMaxFlow::SaturateOutgoingArcsFromSource() { SCOPED_TIME_STAT(&stats_); const NodeIndex num_nodes = graph_->num_nodes(); - // If sink_ or source_ already have kMaxFlowQuantity, then there is no + // If sink_ or source_ already have kMaxFlowSum, then there is no // point pushing more flow since it will cause an integer overflow. - if (node_excess_[sink_] == kMaxFlowQuantity) return false; - if (node_excess_[source_] == -kMaxFlowQuantity) return false; + if (node_excess_[sink_] == kMaxFlowSum) return false; + if (node_excess_[source_] == -kMaxFlowSum) return false; bool flow_pushed = false; for (const ArcIndex arc : graph_->OutgoingArcs(source_)) { - const FlowType flow = residual_arc_capacity_[arc]; + const ArcFlowType flow = residual_arc_capacity_[arc]; // This is a special IsAdmissible() condition for the source. if (flow == 0 || node_potential_[Head(arc)] >= num_nodes) continue; // We are careful in case the sum of the flow out of the source is greater - // than kMaxFlowQuantity to avoid overflow. - const FlowType current_flow_out_of_source = -node_excess_[source_]; + // than kMaxFlowSum to avoid overflow. + const FlowSumType current_flow_out_of_source = -node_excess_[source_]; DCHECK_GE(flow, 0) << flow; DCHECK_GE(current_flow_out_of_source, 0) << current_flow_out_of_source; - const FlowType capped_flow = kMaxFlowQuantity - current_flow_out_of_source; - if (capped_flow < flow) { + const FlowSumType capped_flow = kMaxFlowSum - current_flow_out_of_source; + if (capped_flow < static_cast(flow)) { // We push as much flow as we can so the current flow on the network will - // be kMaxFlowQuantity. + // be kMaxFlowSum. // Since at the beginning of the function, current_flow_out_of_source - // was different from kMaxFlowQuantity, we are sure to have pushed some + // was different from kMaxFlowSum, we are sure to have pushed some // flow before if capped_flow is 0. if (capped_flow == 0) return true; - PushFlow(capped_flow, source_, arc); + PushFlow(static_cast(capped_flow), source_, arc); return true; } PushFlow(flow, source_, arc); @@ -1106,9 +1137,10 @@ bool GenericMaxFlow::SaturateOutgoingArcsFromSource() { return flow_pushed; } -template -void GenericMaxFlow::PushFlow(FlowType flow, NodeIndex tail, - ArcIndex arc) { +template +void GenericMaxFlow::PushFlow(ArcFlowType flow, + NodeIndex tail, + ArcIndex arc) { SCOPED_TIME_STAT(&stats_); // Update the residual capacity of the arc and its opposite arc. @@ -1125,28 +1157,34 @@ void GenericMaxFlow::PushFlow(FlowType flow, NodeIndex tail, // however that we cannot check this because when we cancel the flow on a // cycle in PushFlowExcessBackToSource(), we may break this invariant during // the operation even if it is still valid at the end. - node_excess_[tail] -= flow; - node_excess_[Head(arc)] += flow; + node_excess_[tail] -= static_cast(flow); + node_excess_[Head(arc)] += static_cast(flow); } -template -void GenericMaxFlow::FastPushFlow(FlowType flow, - NodeIndex tail, ArcIndex arc, - FlowType* node_excess) { +template +void GenericMaxFlow::PushAsMuchFlowAsPossible( + NodeIndex tail, ArcIndex arc, FlowSumType* node_excess) { SCOPED_TIME_STAT(&stats_); + // We either push all node_excess or saturate the arc by consuming all its + // residual capacity. + const ArcFlowType flow = static_cast( + std::min(static_cast(residual_arc_capacity_[arc]), + node_excess[tail])); + DCHECK_NE(flow, 0); residual_arc_capacity_[arc] -= flow; residual_arc_capacity_[Opposite(arc)] += flow; DCHECK_GE(residual_arc_capacity_[arc], 0); DCHECK_GE(residual_arc_capacity_[Opposite(arc)], 0); - node_excess[tail] -= flow; - node_excess[Head(arc)] += flow; + node_excess[tail] -= static_cast(flow); + node_excess[Head(arc)] += static_cast(flow); } -template -void GenericMaxFlow::InitializeActiveNodeContainer() { +template +void GenericMaxFlow::InitializeActiveNodeContainer() { SCOPED_TIME_STAT(&stats_); DCHECK(IsEmptyActiveNodeContainer()); const NodeIndex num_nodes = graph_->num_nodes(); @@ -1161,8 +1199,8 @@ void GenericMaxFlow::InitializeActiveNodeContainer() { } } -template -void GenericMaxFlow::RefineWithGlobalUpdate() { +template +void GenericMaxFlow::RefineWithGlobalUpdate() { SCOPED_TIME_STAT(&stats_); const NodeIndex num_nodes = graph_->num_nodes(); @@ -1170,15 +1208,15 @@ void GenericMaxFlow::RefineWithGlobalUpdate() { // Usually SaturateOutgoingArcsFromSource() will saturate all the arcs from // the source in one go, and we will loop just once. But in case we can push - // more than kMaxFlowQuantity out of the source the loop is as follow: - // - Push up to kMaxFlowQuantity out of the source on the admissible outgoing + // more than kMaxFlowSum out of the source the loop is as follow: + // - Push up to kMaxFlowSum out of the source on the admissible outgoing // arcs. Stop if no flow was pushed. // - Compute the current max-flow. This will push some flow back to the // source and render more outgoing arcs from the source not admissible. // // TODO(user): This may not be the most efficient algorithm if we need to loop // many times. An alternative may be to handle the source like the other nodes - // in the algorithm, initially putting an excess of kMaxFlowQuantity on it, + // in the algorithm, initially putting an excess of kMaxFlowSum on it, // and making the source active like any other node with positive excess. To // investigate. while (SaturateOutgoingArcsFromSource()) { @@ -1229,13 +1267,14 @@ void GenericMaxFlow::RefineWithGlobalUpdate() { } } -template -void GenericMaxFlow::Discharge(const NodeIndex node) { +template +void GenericMaxFlow::Discharge( + const NodeIndex node) { SCOPED_TIME_STAT(&stats_); const NodeIndex num_nodes = graph_->num_nodes(); // We cache this as this is a hot-loop and we don't want bound checking. - FlowType* const node_excess = node_excess_.get(); + FlowSumType* const node_excess = node_excess_.get(); NodeHeight* const node_potentials = node_potential_.get(); ArcIndex* const first_admissible_arc = first_admissible_arc_.get(); @@ -1252,9 +1291,7 @@ void GenericMaxFlow::Discharge(const NodeIndex node) { // push the sink_, but that is handled properly in Refine(). PushActiveNode(head); } - const FlowType delta = - std::min(node_excess[node], residual_arc_capacity_[arc]); - FastPushFlow(delta, node, arc, node_excess); + PushAsMuchFlowAsPossible(node, arc, node_excess); if (node_excess[node] == 0) { first_admissible_arc[node] = arc; // arc may still be admissible. return; @@ -1269,8 +1306,8 @@ void GenericMaxFlow::Discharge(const NodeIndex node) { } } -template -void GenericMaxFlow::Relabel(NodeIndex node) { +template +void GenericMaxFlow::Relabel(NodeIndex node) { SCOPED_TIME_STAT(&stats_); // Because we use a relaxed version, this is no longer true if the // first_admissible_arc_[node] was not actually the first arc! @@ -1300,25 +1337,26 @@ void GenericMaxFlow::Relabel(NodeIndex node) { first_admissible_arc_[node] = first_admissible_arc; } -template -typename Graph::ArcIndex GenericMaxFlow::Opposite( +template +typename Graph::ArcIndex GenericMaxFlow::Opposite( ArcIndex arc) const { return graph_->OppositeArc(arc); } -template -bool GenericMaxFlow::IsArcDirect(ArcIndex arc) const { +template +bool GenericMaxFlow::IsArcDirect( + ArcIndex arc) const { return IsArcValid(arc) && arc >= 0; } -template -bool GenericMaxFlow::IsArcValid(ArcIndex arc) const { +template +bool GenericMaxFlow::IsArcValid(ArcIndex arc) const { return graph_->IsArcValid(arc); } -template +template template -void GenericMaxFlow::ComputeReachableNodes( +void GenericMaxFlow::ComputeReachableNodes( NodeIndex start, std::vector* result) { // If start is not a valid node index, it can reach only itself. // Note(user): This is needed because source and sink are given independently @@ -1349,8 +1387,8 @@ void GenericMaxFlow::ComputeReachableNodes( *result = bfs_queue_; } -template -FlowModelProto GenericMaxFlow::CreateFlowModel() { +template +FlowModelProto GenericMaxFlow::CreateFlowModel() { FlowModelProto model; model.set_problem_type(FlowModelProto::MAX_FLOW); for (int n = 0; n < graph_->num_nodes(); ++n) { diff --git a/ortools/graph/generic_max_flow_test.cc b/ortools/graph/generic_max_flow_test.cc index 99d22fead8..20e82e6e7f 100644 --- a/ortools/graph/generic_max_flow_test.cc +++ b/ortools/graph/generic_max_flow_test.cc @@ -19,6 +19,7 @@ #include #include #include +#include #include #include #include @@ -69,7 +70,7 @@ typename GenericMaxFlow::Status MaxFlowTester( } EXPECT_TRUE(max_flow.Solve()); if (max_flow.status() == GenericMaxFlow::OPTIMAL) { - const typename GenericMaxFlow::FlowQuantityT total_flow = + const typename GenericMaxFlow::FlowSumType total_flow = max_flow.GetOptimalFlow(); EXPECT_EQ(expected_total_flow, total_flow); for (int arc = 0; arc < num_arcs; ++arc) { @@ -177,7 +178,7 @@ TYPED_TEST(GenericMaxFlowTest, HugeCapacity) { } TYPED_TEST(GenericMaxFlowTest, FlowQuantityOverflowLimitCase) { - using FlowQuantity = typename GenericMaxFlow::FlowQuantityT; + using FlowQuantity = typename GenericMaxFlow::FlowSumType; const FlowQuantity kCapacityMax = std::numeric_limits::max(); const FlowQuantity kHalfLow = kCapacityMax / 2; const FlowQuantity kHalfHigh = kCapacityMax - kHalfLow; @@ -199,7 +200,7 @@ TYPED_TEST(GenericMaxFlowTest, FlowQuantityOverflowLimitCase) { } TYPED_TEST(GenericMaxFlowTest, FlowQuantityOverflow) { - using FlowQuantity = typename GenericMaxFlow::FlowQuantityT; + using FlowQuantity = typename GenericMaxFlow::FlowSumType; const FlowQuantity kCapacityMax = std::numeric_limits::max(); const int kNumNodes = 4; const int kNumArcs = 4; @@ -264,6 +265,79 @@ TYPED_TEST(GenericMaxFlowTest, FlowOnDisconnectedGraph2) { kExpectedTotalFlow, &source_cut, &sink_cut)); } +// Using a custom type should allow to verify that we didn't forget a conversion +// somewhere. +// +// TODO(user): Unfortunately, there is no open-source version of strong int +// supporting uint16_t... +struct StrongUint16 { + constexpr StrongUint16() : v(0) {} + constexpr StrongUint16(int v) : v(v) {} + + explicit operator int64_t() const { return static_cast(v); } + + StrongUint16 operator-() const { return StrongUint16(-v); } + void operator+=(StrongUint16 o) { v += o.v; } + void operator-=(StrongUint16 o) { v -= o.v; } + + uint16_t v; +}; +std::ostream& operator<<(std::ostream& os, StrongUint16 a) { return os << a.v; } +constexpr StrongUint16 operator+(StrongUint16 a, StrongUint16 b) { + return StrongUint16(a.v + b.v); +} +constexpr StrongUint16 operator-(StrongUint16 a, StrongUint16 b) { + return StrongUint16(a.v - b.v); +} +constexpr bool operator<=(StrongUint16 a, StrongUint16 b) { return a.v <= b.v; } +constexpr bool operator>=(StrongUint16 a, StrongUint16 b) { return a.v >= b.v; } +constexpr bool operator>(StrongUint16 a, StrongUint16 b) { return a.v > b.v; } +constexpr bool operator<(StrongUint16 a, StrongUint16 b) { return a.v < b.v; } +constexpr bool operator==(StrongUint16 a, StrongUint16 b) { return a.v == b.v; } + +// This is a bit hacky, but we need to define the overload in the std namespace. +} // namespace +} // namespace operations_research +namespace std { +template <> +struct numeric_limits { + static constexpr operations_research::StrongUint16 max() { + return operations_research::StrongUint16(numeric_limits::max()); + } +}; +} // namespace std +namespace operations_research { +namespace { + +TYPED_TEST(GenericMaxFlowTest, SmallFlowTypes) { + absl::BitGen random; + const int num_nodes = 1'000; + const int num_arcs = num_nodes * num_nodes; + + // Lets generate and build a random graph. + // We should have more than enough arc to make it fully connected. + TypeParam graph(num_nodes, num_arcs); + for (int arc = 0; arc < num_arcs; ++arc) { + graph.AddArc(absl::Uniform(random, 0, num_nodes), + absl::Uniform(random, 0, num_nodes)); + } + graph.Build(); + + typedef GenericMaxFlow MaxFlowA; + typedef GenericMaxFlow MaxFlowB; + MaxFlowA max_flow_A(&graph, /*source=*/0, /*sink=*/num_nodes - 1); + MaxFlowB max_flow_B(&graph, /*source=*/0, /*sink=*/num_nodes - 1); + for (int arc = 0; arc < num_arcs; ++arc) { + const uint16_t capa = + absl::Uniform(random, 0, std::numeric_limits::max() / 2); + max_flow_A.SetArcCapacity(arc, capa); + max_flow_B.SetArcCapacity(arc, capa); + } + EXPECT_EQ(max_flow_A.Solve(), MaxFlowA::OPTIMAL); + EXPECT_EQ(max_flow_B.Solve(), MaxFlowB::OPTIMAL); + EXPECT_EQ(max_flow_A.GetOptimalFlow(), max_flow_B.GetOptimalFlow()); +} + template void AddSourceAndSink(const typename Graph::NodeIndex num_tails, const typename Graph::NodeIndex num_heads, Graph* graph) { diff --git a/ortools/graph/shortest_paths.h b/ortools/graph/shortest_paths.h index 80a964bf60..07e4c9908f 100644 --- a/ortools/graph/shortest_paths.h +++ b/ortools/graph/shortest_paths.h @@ -184,14 +184,6 @@ class GenericPathContainer { std::unique_ptr container_; }; -// TODO(b/385094969): Remove this alias when all clients are migrated. -class LegacyIgnoredGraphType { - public: - using NodeIndex = int32_t; - static constexpr NodeIndex kNilNode = std::numeric_limits::max(); -}; -using PathContainer = GenericPathContainer; - // Utility function which returns a vector containing all nodes of a graph. template void GetGraphNodes(const GraphType& graph,