121#ifndef OR_TOOLS_GRAPH_ONE_TREE_LOWER_BOUND_H_
122#define OR_TOOLS_GRAPH_ONE_TREE_LOWER_BOUND_H_
133#include "absl/types/span.h"
157template <
typename CostType>
158class VolgenantJonkerEvaluator {
161 : step1_initialized_(false),
165 : MaxIterations(number_of_nodes)),
166 number_of_nodes_(number_of_nodes) {}
168 bool Next() {
return iteration_++ < max_iterations_; }
171 return (1.0 * (iteration_ - 1) * (2 * max_iterations_ - 5) /
172 (2 * (max_iterations_ - 1))) *
174 (iteration_ - 2) * step1_ +
175 (0.5 * (iteration_ - 1) * (iteration_ - 2) /
176 ((max_iterations_ - 1) * (max_iterations_ - 2))) *
181 absl::Span<const int>
degrees) {
182 if (!step1_initialized_) {
183 step1_initialized_ =
true;
184 UpdateStep(one_tree_cost);
188 void OnNewWMax(CostType one_tree_cost) { UpdateStep(one_tree_cost); }
193 static int MaxIterations(
int number_of_nodes) {
194 return static_cast<int>(28 * std::pow(number_of_nodes, 0.62));
197 void UpdateStep(CostType one_tree_cost) {
198 step1_ = one_tree_cost / (2 * number_of_nodes_);
201 bool step1_initialized_;
204 const int max_iterations_;
205 const int number_of_nodes_;
210template <
typename CostType,
typename CostFunction>
211class HeldWolfeCrowderEvaluator {
213 HeldWolfeCrowderEvaluator(
int number_of_nodes,
const CostFunction& cost)
215 number_of_iterations_(2 * number_of_nodes),
221 ChristofidesPathSolver<CostType, int64_t, int, CostFunction> solver(
222 number_of_nodes, cost);
223 upper_bound_ = solver.TravelingSalesmanCost();
227 const int min_iterations = 2;
228 if (iteration_ >= number_of_iterations_) {
229 number_of_iterations_ /= 2;
230 if (number_of_iterations_ < min_iterations)
return false;
239 double GetStep()
const {
return step_; }
241 void OnOneTree(CostType one_tree_cost,
double w,
242 absl::Span<const int>
degrees) {
245 const double delta = degree - 2;
248 step_ = lambda_ * (upper_bound_ -
w) / norm;
251 void OnNewWMax(CostType one_tree_cost) {}
255 int number_of_iterations_;
256 CostType upper_bound_;
266template <
typename CostFunction>
267std::set<std::pair<int, int>> NearestNeighbors(
int number_of_nodes,
268 int number_of_neighbors,
269 const CostFunction& cost) {
270 using CostType =
decltype(cost(0, 0));
271 std::set<std::pair<int, int>> nearest;
272 for (
int i = 0;
i < number_of_nodes; ++
i) {
273 std::vector<std::pair<CostType, int>> neighbors;
274 neighbors.reserve(number_of_nodes - 1);
275 for (
int j = 0; j < number_of_nodes; ++j) {
277 neighbors.emplace_back(cost(i, j), j);
280 int size = neighbors.size();
281 if (number_of_neighbors < size) {
282 std::nth_element(neighbors.begin(),
283 neighbors.begin() + number_of_neighbors - 1,
285 size = number_of_neighbors;
287 for (
int j = 0; j < size; ++j) {
288 nearest.insert({
i, neighbors[j].second});
289 nearest.insert({neighbors[j].second,
i});
297template <
typename CostFunction>
298void AddArcsFromMinimumSpanningTree(
int number_of_nodes,
299 const CostFunction& cost,
300 std::set<std::pair<int, int>>* arcs) {
302 const std::vector<int> mst =
304 return cost(graph.Tail(
arc), graph.Head(
arc));
306 for (
int arc : mst) {
307 arcs->insert({graph.Tail(
arc), graph.Head(
arc)});
308 arcs->insert({graph.Head(
arc), graph.Tail(
arc)});
314template <
typename CostFunction,
typename GraphType,
typename AcceptFunction>
315int GetNodeMinimizingEdgeCostToSource(
const GraphType& graph,
int source,
316 const CostFunction& cost,
317 AcceptFunction accept) {
319 double best_edge_cost = 0;
320 for (
const auto node : graph.AllNodes()) {
322 const double edge_cost = cost(node, source);
323 if (best_node == -1 || edge_cost < best_edge_cost) {
325 best_edge_cost = edge_cost;
335template <
typename CostFunction,
typename GraphType,
typename CostType>
336std::vector<int> ComputeOneTree(
const GraphType& graph,
337 const CostFunction& cost,
338 const std::vector<double>& weights,
339 const std::vector<int>& sorted_arcs,
340 CostType* one_tree_cost) {
341 const auto weighed_cost = [&cost, &weights](
int from,
int to) {
342 return cost(from,
to) + weights[from] + weights[
to];
345 std::vector<int> mst;
346 if (!sorted_arcs.empty()) {
347 mst = BuildKruskalMinimumSpanningTreeFromSortedArcs<GraphType>(graph,
350 mst = BuildPrimMinimumSpanningTree<GraphType>(
351 graph, [&weighed_cost, &graph](
int arc) {
352 return weighed_cost(graph.Tail(
arc), graph.Head(
arc));
355 std::vector<int>
degrees(graph.num_nodes() + 1, 0);
357 for (
int arc : mst) {
360 *one_tree_cost += cost(graph.Tail(
arc), graph.Head(
arc));
364 const int extra_node = graph.num_nodes();
365 const auto update_one_tree = [extra_node, one_tree_cost, &
degrees,
367 *one_tree_cost += cost(node, extra_node);
371 const int node = GetNodeMinimizingEdgeCostToSource(
372 graph, extra_node, weighed_cost,
373 [extra_node](
int n) {
return n != extra_node; });
374 update_one_tree(node);
375 update_one_tree(GetNodeMinimizingEdgeCostToSource(
376 graph, extra_node, weighed_cost,
377 [extra_node, node](
int n) {
return n != extra_node && n != node; }));
382template <
typename CostFunction,
typename Algorithm>
383double ComputeOneTreeLowerBoundWithAlgorithm(
int number_of_nodes,
384 int nearest_neighbors,
385 const CostFunction& cost,
386 Algorithm* algorithm) {
387 if (number_of_nodes < 2)
return 0;
388 if (number_of_nodes == 2)
return cost(0, 1) + cost(1, 0);
389 using CostType =
decltype(cost(0, 0));
390 auto nearest = NearestNeighbors(number_of_nodes - 1, nearest_neighbors, cost);
394 AddArcsFromMinimumSpanningTree(number_of_nodes - 1, cost, &nearest);
396 for (
const auto&
arc : nearest) {
397 graph.AddArc(
arc.first,
arc.second);
399 std::vector<double> weights(number_of_nodes, 0);
400 std::vector<double> best_weights(number_of_nodes, 0);
401 double max_w = -std::numeric_limits<double>::infinity();
404 while (algorithm->Next()) {
405 CostType one_tree_cost = 0;
406 const std::vector<int>
degrees =
407 ComputeOneTree(graph, cost, weights, {}, &one_tree_cost);
408 algorithm->OnOneTree(one_tree_cost,
w,
degrees);
410 for (
int j = 0; j < number_of_nodes; ++j) {
415 best_weights = weights;
416 algorithm->OnNewWMax(one_tree_cost);
418 const double step = algorithm->GetStep();
419 for (
int j = 0; j < number_of_nodes; ++j) {
420 weights[j] += step * (
degrees[j] - 2);
427 CostType one_tree_cost = 0;
431 const std::vector<int>
degrees =
432 ComputeOneTree(complete_graph, cost, best_weights, {}, &one_tree_cost);
434 for (
int j = 0; j < number_of_nodes; ++j) {
435 w += best_weights[j] * (
degrees[j] - 2);
441struct TravelingSalesmanLowerBoundParameters {
447 Algorithm algorithm = VolgenantJonker;
450 int volgenant_jonker_iterations = 0;
452 int nearest_neighbors = 40;
456template <
typename CostFunction>
457double ComputeOneTreeLowerBoundWithParameters(
458 int number_of_nodes,
const CostFunction& cost,
459 const TravelingSalesmanLowerBoundParameters&
parameters) {
460 using CostType =
decltype(cost(0, 0));
462 case TravelingSalesmanLowerBoundParameters::VolgenantJonker: {
463 VolgenantJonkerEvaluator<CostType> algorithm(
464 number_of_nodes,
parameters.volgenant_jonker_iterations);
465 return ComputeOneTreeLowerBoundWithAlgorithm(
466 number_of_nodes,
parameters.nearest_neighbors, cost, &algorithm);
469 case TravelingSalesmanLowerBoundParameters::HeldWolfeCrowder: {
470 HeldWolfeCrowderEvaluator<CostType, CostFunction> algorithm(
471 number_of_nodes, cost);
472 return ComputeOneTreeLowerBoundWithAlgorithm(
473 number_of_nodes,
parameters.nearest_neighbors, cost, &algorithm);
476 LOG(ERROR) <<
"Unsupported algorithm: " <<
parameters.algorithm;
484template <
typename CostFunction>
485double ComputeOneTreeLowerBound(
int number_of_nodes,
const CostFunction& cost) {
486 TravelingSalesmanLowerBoundParameters
parameters;
487 return ComputeOneTreeLowerBoundWithParameters(number_of_nodes, cost,
In SWIG mode, we don't want anything besides these top-level includes.
std::vector< typename Graph::ArcIndex > BuildPrimMinimumSpanningTree(const Graph &graph, const ArcValue &arc_value)
void OnOneTree(CostType one_tree_cost, double w, absl::Span< const int > degrees)
trees with all degrees equal w the current value of w
void OnNewWMax(CostType one_tree_cost)
trees with all degrees equal to
trees with all degrees equal w the current value of int max_iterations
trees with all degrees equal w the current value of degrees