最小カット問題のk値への一般化 (src/Optimization/monge_mincut.hpp)

View this file on GitHub
View document part on GitHub
Last update: 2023-01-21 23:17:22+09:00
Include: #include "src/Optimization/monge_mincut.hpp"

phi関数等で ∞ を返すときはその大きさに注意 (大きすぎるとオーバーフロー)

関数

名前	概要
`monge_mincut<MF>(n,k,theta,phi)`	$\begin{aligned} \min_{x\in\lbrace 0,\dots,k-1 \rbrace^n} \left(\sum_{i=0}^{n-1} \theta_i(x_i) + \sum_{i\lt j}\phi_{ij}(x_i,x_j)\right) \end{aligned}$ を求める(その時の解$x$も求める). ただし $\phi_{ij}(\cdot,\cdot)$ はMonge ( 特に $\phi_{ij}(a,b)+\phi_{ij}(a+1,b+1)\leq \phi_{ij}(a,b+1)+\phi_{ij}(a+1,b)$ ) テンプレート引数で最大流クラスを渡す.

使用例

 auto theta= [&](int i, int xi) -> long long {
  ...
 };
 auto phi= [&](int i, int j, int xi, int xj) -> long long {
  ...
 };
 using MF = MaxFlow<Dinic<long long>>;
 auto [ans, x] = monge_mincut<MF>(N, k, theta, phi);

Verify

参考

https://noshi91.hatenablog.com/entry/2021/06/29/044225

Verified with

Code

#pragma once
#include <vector>
#include <limits>
#include <cassert>
#include <utility>
template <typename MF, typename Th, typename Ph> auto monge_mincut(int n, int k, Th theta, Ph phi) {
 using T= typename MF::flow_t;
 static constexpr T INF= std::numeric_limits<T>::max();
 T ret= 0;
 MF graph;
 int s= graph.add_vertex(), t= graph.add_vertex();
 std::vector<std::vector<int>> x;
 std::vector<std::vector<T>> th(n, std::vector<T>(k));
 for (int i= 0; i < n; i++) {
  x.emplace_back(graph.add_vertices(k - 1));
  for (int l= 1; l < k - 1; l++) graph.add_edge(x[i][l], x[i][l - 1], INF);
  for (int l= 0; l < k; l++) th[i][l]= theta(i, l);
 }
 for (int i= 0; i < n; i++)
  for (int j= i + 1; j < n; j++) {
   for (int l= 0; l < k - 1; l++)
    for (int m= 0; m < k - 1; m++) {
     T cost= -phi(i, j, l + 1, m + 1) + phi(i, j, l, m + 1) + phi(i, j, l + 1, m) - phi(i, j, l, m);
     assert(cost >= 0);  // monge
     if (cost > 0) graph.add_edge(x[i][l], x[j][m], cost);
    }
   for (int l= 0; l < k; l++) th[i][l]+= phi(i, j, l, k - 1);
   for (int l= 0; l < k; l++) th[j][l]+= phi(i, j, 0, l);
   ret-= phi(i, j, 0, k - 1);
  }
 for (int i= 0; i < n; i++) {
  ret+= th[i][0];
  for (int l= 0; l < k - 1; l++) {
   T cost= th[i][l] - th[i][l + 1];
   if (cost > 0) graph.add_edge(s, x[i][l], cost), ret-= cost;
   if (cost < 0) graph.add_edge(x[i][l], t, -cost);
  }
 }
 ret+= graph.maxflow(s, t);
 auto y= graph.mincut(s);
 std::vector<int> sol(n, k - 1);
 for (int i= 0; i < n; i++)
  for (int l= 0; l < k - 1; l++)
   if (!y[x[i][l]]) sol[i]= l, l= k;
 return std::make_pair(ret, sol);
}

#line 2 "src/Optimization/monge_mincut.hpp"
#include <vector>
#include <limits>
#include <cassert>
#include <utility>
template <typename MF, typename Th, typename Ph> auto monge_mincut(int n, int k, Th theta, Ph phi) {
 using T= typename MF::flow_t;
 static constexpr T INF= std::numeric_limits<T>::max();
 T ret= 0;
 MF graph;
 int s= graph.add_vertex(), t= graph.add_vertex();
 std::vector<std::vector<int>> x;
 std::vector<std::vector<T>> th(n, std::vector<T>(k));
 for (int i= 0; i < n; i++) {
  x.emplace_back(graph.add_vertices(k - 1));
  for (int l= 1; l < k - 1; l++) graph.add_edge(x[i][l], x[i][l - 1], INF);
  for (int l= 0; l < k; l++) th[i][l]= theta(i, l);
 }
 for (int i= 0; i < n; i++)
  for (int j= i + 1; j < n; j++) {
   for (int l= 0; l < k - 1; l++)
    for (int m= 0; m < k - 1; m++) {
     T cost= -phi(i, j, l + 1, m + 1) + phi(i, j, l, m + 1) + phi(i, j, l + 1, m) - phi(i, j, l, m);
     assert(cost >= 0);  // monge
     if (cost > 0) graph.add_edge(x[i][l], x[j][m], cost);
    }
   for (int l= 0; l < k; l++) th[i][l]+= phi(i, j, l, k - 1);
   for (int l= 0; l < k; l++) th[j][l]+= phi(i, j, 0, l);
   ret-= phi(i, j, 0, k - 1);
  }
 for (int i= 0; i < n; i++) {
  ret+= th[i][0];
  for (int l= 0; l < k - 1; l++) {
   T cost= th[i][l] - th[i][l + 1];
   if (cost > 0) graph.add_edge(s, x[i][l], cost), ret-= cost;
   if (cost < 0) graph.add_edge(x[i][l], t, -cost);
  }
 }
 ret+= graph.maxflow(s, t);
 auto y= graph.mincut(s);
 std::vector<int> sol(n, k - 1);
 for (int i= 0; i < n; i++)
  for (int l= 0; l < k - 1; l++)
   if (!y[x[i][l]]) sol[i]= l, l= k;
 return std::make_pair(ret, sol);
}