Hoi_koro memo

Editorialと異なる解法をしたら更新するblog

AtCoder 全国統一プログラミング王決定戦本戦/NIKKEI Programming Contest 2019 H - Homework Scheduling

atcoder.jp

問題概要

$ N $個の宿題があり、1日に1つずつ終わらせていく。各宿題には期限$ A_i $ が設定されており、宿題$ i $をDay $ A_i $までに終わらせると$ X_i $点、それ以後に終わらせると$ Y_i $点が得られる。$ 1\leq k \leq N $なる$ k $について、Day 1 ~ Day $ k $で得られる得点を最大化せよ。

制約

  • $1\leq N\leq 2\times 10^5$
  • $1 \leq A_i \leq N $
  • $ 1\leq Y_i < X_i \leq N $

解法

はじめに$A_i$が昇順になるように宿題を並べ替える。

最小費用流への言い換え

グラフの構築

この問題のDay $k$に関する答えは、次のグラフにおける流量$ k $の$S\to T$最小費用流のコストの-1倍である。

f:id:hiko3r:20190227225828p:plain:w600

ただし、$ Q_i \to T$の辺たちは次の方法によってはられる。

「$i = 1,2,\dots,N$のそれぞれに対して、$A_k< i$となる$ k $の個数を$c_i$とする。$Q_{c_i+1} \to T$に容量1の辺をはる。」

$ Q_i \to T$の辺に左から順に$r_1,\dots,r_n$と名前をつけておく。

グラフの例

問題の例2について、$\{A_i\}_{i=1}^6=\{1,2,2,3,4,4\}$なので$\{c_i\} _{i=1}^6= \{1,3,4,6,7,7\}$となる。 グラフは下のようになる。

f:id:hiko3r:20190228003004p:plain:w480

グラフの最小費用流へ帰着できることの説明

グラフ上での$ S\to T$パスは、

  • {S\to P_i \to Q_i \to Q_j \overset{r_k}{\rightarrow} T} ( $j \leq i$)
  • $S\to P_i \to Q_{N+1} \to Q_j \overset{r_k}{\rightarrow} T$ ( $j\leq N+1$)

のどちらかの形になる。

前者を宿題$i$を$k$日目に終わらせ、かつ$X_i$点を獲得した(すなわち$k\leq A_i$)ことに対応させる。また、後者を宿題$i$を$k$日目に終わらせ、かつ$Y_i$点を獲得した(すなわち$k> A_i$)ことに対応させる。 従って、下図青のパスは宿題2を1日目に終わらせ、$X_2$点を得たことに対応し、赤のパスは宿題3を4日目に終わらせ、$Y_3$点を得たこと に対応している。

f:id:hiko3r:20190228003009p:plain:w480

ここで、

(A)任意のDay $k$までの宿題の終わらせ方がグラフ上の容量$k$のフローで表現でき、かつ( 宿題の得点の和 )=( フローのコスト )$\times (-1) $であること

逆に(B)グラフ上の任意の容量$k$のフローがDay $k$までの宿題の終わらせ方を表し、かつ( 宿題の得点の和 )=( フローのコスト )$\times (-1) $であること

がわかれば、最小費用流へ問題を言い換えてよいと言える。それぞれが正しいことは次のように考えるとわかる。

(A) Day $k$までのある宿題の終わらせ方が与えられたとする。日ごとに終わらせた宿題をグラフ上の$S\to T$パスに置き換えて、$k$本の$S\to T$パスを作る。$r_1,\dots,r_n$の作り方から、$k$本のパスに流量1を流したフローがグラフの容量制約をみたすことがわかる。

(B) グラフ上の容量$k$のフローが与えられたとする。辺$r_1,\dots,r_n$が左から優先して使われるようにパスを左に詰める。たとえば上図の容量2のフローの場合、赤いパスを{S\to P_3 \to Q_7 \to Q_2 \overset{r_2}{\rightarrow} T}にずらす。この変更はフローのコストを変えない。変更を行った後でそれぞれのパスを宿題の終わらせ方に置き換える。辺$r_1,\dots,r_n$の作り方から、制約をみたす宿題の終わらせ方であることがわかる。

最小費用流問題の解き方

ad-hocに工夫したPrimal-Dualによって解く。 (Primal-Dualの参考 Spaghetti Source - 最小費用流 (Primal-Dual) 最小費用流(Primal-Dual) | Luzhiled’s memo

Primal-Dualは、残余グラフ上の最小コストの$S\to T$パスに流量1のフローを流す操作を繰り返していく方法である。この問題ではグラフの構造から最小コストの$S\to T$パスを求める部分が高速にできるのがポイントである。

なお、editorial中にある仮定

  • 辺 Qi → T の容量が正であるすべての i に対して、辺 Qi−1 → Qi の容量が 0 である。

というのは、前項とは逆に辺$r_1,\dots,r_n$が右から優先して使われるようにパスを右に詰めることを意味している。この記事中でもそのように辺$r_1,\dots,r_n$を右から使うことにする。

残余グラフ上の$S\to T$パスは次の4種類に分けられる。

  • A : {S\to P_i \to Q_i \to Q_j \overset{r_l}{\rightarrow} T} ( $j \leq i$)
  • B : $S\to P_i \to Q_{N+1} \to Q_j \overset{r_l}{\rightarrow} T$ ( $j\leq N+1$)
  • C : $S\to P_i \to Q_i \to Q_j \to P_j \to Q_{N+1} \to Q_k \overset{r_l}{\rightarrow} T$ ( $i> j$)
  • D : $S\to P_i \to Q_i \to Q_j \to P_j \to Q_{N+1} \to Q_k \overset{r_l}{\rightarrow} T$ ( $i< j$)

仮定の部分でも述べたように、$r_l$は残余グラフ上で空いている中でできるだけ右、つまり$l$の大きいものを使う。従って、A~Dのタイプの中での最小コストは$i,j$を決めることで一意に定まる。

C, Dタイプの最小コストを与える$i,j$は、解説にあるSegment Treeを実装すると求められる。Cタイプの最小コストは$C_i+D_j ~(i>j)$の最小値である。Dタイプの最小コストは$m\neq 0$ ならば単に$C_i+D_j ~(i<j)$の最小値であり、$m=0$ならば$C_i+D_j ~(i<j, E_i,E_j >m)$の最小値である。

Aタイプについては,点更新と区間minクエリに対応したSegmentTreeを用意して、 {S\to P_i \to Q_i \to Q_j \overset{r_l}{\rightarrow} T} のパスがあるなら$seg[i]=$そのコスト、ないなら$seg[i]=\infty$となるように更新していけば求められる。Bタイプも同様。

実装の際には、辺$r_1,\dots,r_n$を右から順に使うという制約をみたすため、残余グラフ上の$ Q_i \to T$の辺のうち容量が正のものを管理しておく必要がある。

アルゴリズム

$1\leq k\leq N$について,

  • A,B,C,Dそれぞれのタイプでの最小コストの$S\to T$パスを求める
  • その中で最小コストの$S\to T$パスに流量1を流し、残余グラフを更新する
  • 現時点でのコストの(-1)倍を答える

を繰り返す。

解答例

#include <bits/stdc++.h>

#define DBG(...) ;
using LL = long long;
constexpr LL LINF = 334ll << 50;
constexpr int INF = 15 << 26;

namespace Problem {
using namespace std;
struct Min_Id {
  LL val;
  int id;
  Min_Id() : val(LINF), id(-1){};
  Min_Id(LL v, const int id) : val(v), id(id){};
  bool operator<(const Min_Id &b) const {
    return val < b.val;
  }
  bool operator==(const Min_Id &b) const {
    return val == b.val;
  }
};

struct Min_IdLR {
  LL val;
  int idl, idr;
  Min_IdLR() : val(LINF * 2), idl(INF), idr(-1){};
  Min_IdLR(LL v, const int l, const int r) : val(v), idl(l), idr(r){};
  bool operator<(const Min_IdLR &b) const {
    return val < b.val;
  }
  bool operator==(const Min_IdLR &b) const {
    return val == b.val;
  }
  tuple<const LL &, const int &, const int &> ref() const {
    return tie(val, idl, idr);
  }
};

struct Min_Range {
  LL val;
  int idl, idr;
  Min_Range() : val(LINF), idl(INF), idr(-1){};
  Min_Range(LL v, const int l, const int r) : val(v), idl(l), idr(r){};
  bool operator<(const Min_Range &b) const {
    return val < b.val;
  }
  bool operator==(const Min_Range &b) const {
    return val == b.val;
  }
};

struct Node {
  Min_Range E;
  Min_Id CL, CR, DL, DR;
  Min_IdLR Dc, CdMinE, CdNoMinE;
  Node() : E(), CL(), CR(), DL(), DR(), Dc(), CdMinE(), CdNoMinE(){};
};

struct Operator {
  LL addC, addD;
  int addE;
  Operator() : addC(0), addD(0), addE(0){};
  Operator(LL c, LL d, int e) : addC(c), addD(d), addE(e){};
};

struct RangeNodeMin {
  using T = Node;
  T operator()(const T &a, const T &b) {
    T ret = T();
    if (a.E < b.E) {
      ret.E = a.E;
      ret.CL = a.CL;
      ret.CR = min({a.CR, b.CL, b.CR});
      ret.DL = a.DL;
      ret.DR = min({a.DR, b.DL, b.DR});
      auto CdMinETrans = Min_IdLR(a.CL.val + min(b.DL, b.DR).val, a.CL.id, min(b.DL, b.DR).id);
      ret.CdMinE = min(a.CdMinE, CdMinETrans);
      auto CdNoMinETrans = Min_IdLR(a.CR.val + min(b.DL, b.DR).val, a.CR.id, min(b.DL, b.DR).id);
      ret.CdNoMinE = min({a.CdNoMinE, b.CdNoMinE, b.CdMinE, CdNoMinETrans});
    } else if (b.E < a.E) {
      ret.E = b.E;
      ret.CL = min({a.CL, a.CR, b.CL});
      ret.CR = b.CR;
      ret.DL = min({a.DL, a.DR, b.DL});
      ret.DR = b.DR;
      auto CdMinETrans = Min_IdLR(b.DR.val + min(a.CL, a.CR).val, min(a.CL, a.CR).id, b.DR.id);
      ret.CdMinE = min(b.CdMinE, CdMinETrans);
      auto CdNoMinETrans = Min_IdLR(b.DL.val + min(a.CL, a.CR).val, min(a.CL, a.CR).id, b.DL.id);
      ret.CdNoMinE = min({a.CdNoMinE, a.CdMinE, b.CdNoMinE, CdNoMinETrans});
    } else {
      ret.E = a.E;
      ret.E.idr = b.E.idr;
      ret.CL = min({a.CL, a.CR, b.CL});
      ret.CR = b.CR;
      ret.DL = a.DL;
      ret.DR = min({a.DR, b.DL, b.DR});
      auto CdMinETrans = Min_IdLR(a.CL.val + min(b.DL, b.DR).val, a.CL.id, min(b.DL, b.DR).id);
      auto CdMinETrans2 = Min_IdLR(b.DR.val + min(a.CL, a.CR).val, min(a.CL, a.CR).id, b.DR.id);
      ret.CdMinE = min({a.CdMinE, b.CdMinE, CdMinETrans, CdMinETrans2});
      auto CdNoMinETrans = Min_IdLR(a.CR.val + b.DL.val, a.CR.id, b.DL.id);
      ret.CdNoMinE = min({a.CdNoMinE, b.CdNoMinE, CdNoMinETrans});
    }
    Min_IdLR minDcTrans = Min_IdLR(min(a.DL.val, a.DR.val) + min(b.CL.val, b.CR.val), min(a.DL, a.DR).id, min(b.CL, b.CR).id);
    ret.Dc = min({a.Dc, b.Dc, minDcTrans});
    return ret;
  }
  static T id() {
    return T();
  };
};

struct RangeAdd {
  using T = Operator;
  T operator()(const T &a, const T &b) {
    Operator ret = a;
    ret.addC += b.addC;
    ret.addD += b.addD;
    ret.addE += b.addE;
    return ret;
  };
  static T id() { return T(); };
};

struct Min_Add {
  using V = RangeNodeMin;
  using A = RangeAdd;
  V::T operator()(const V::T &a, const A::T &b, const int &h) {
    V::T ret = a;
    ret.CL.val += b.addC;
    ret.CR.val += b.addC;
    ret.DL.val += b.addD;
    ret.DR.val += b.addD;
    ret.E.val += b.addE;
    return ret;
  };
};

struct RangeMin {
  using T = pair<long long, int>;
  T operator()(const T &a, const T &b) { return min(a, b); }
  static constexpr T id() { return T(LINF, -1); };
};

struct RangeSet {
  using T = long long;
  T operator()(const T &a, const T &b) { return b; };
  static constexpr T id() { return T(LINF * 2); };
};

struct Min_Set {
  using V = RangeMin;
  using A = RangeSet;
  V::T operator()(const V::T &a, const A::T &b, const int &h) {
    if (b == LINF * 2)
      return a;
    else
      return V::T(b, a.second);
  };
};

template <typename M>
struct LazySegmentTree {
  using Val = typename M::V::T;
  using Del = typename M::A::T;
  typename M::V calc;
  typename M::A composite;
  M act;
  vector<Val> val;
  vector<Del> delay;

  int n, height;
  explicit LazySegmentTree(int size) : n(1), height(0) {
    while (n < size) {
      n <<= 1;
      ++height;
    }
    val = vector<Val>(n << 1, calc.id());
    delay = vector<Del>(n << 1, composite.id());
  }
  explicit LazySegmentTree(const vector<Val> &v) : n(1), height(0) {
    while (n < v.size()) {
      n <<= 1;
      ++height;
    }
    val = vector<Val>(n << 1, calc.id());
    delay = vector<Del>(n << 1, composite.id());
    copy(v.begin(), v.end(), val.begin() + n);
    for (int i = n - 1; i > 0; --i) {
      val[i] = calc(val[i * 2], val[i * 2 + 1]);
    }
  }
  void propagate(int k, int h) {
    val[k] = act(val[k], delay[k], h);
    if (h > 0) {
      delay[2 * k] = composite(delay[2 * k], delay[k]);
      delay[2 * k + 1] = composite(delay[2 * k + 1], delay[k]);
    }
    delay[k] = composite.id();
  }
  inline void update(int l, int r, Del v) { update(l, r, v, 1, 0, n, height); }
  void update(int l, int r, Del v, int k, int a, int b, int h) {
    if (l <= a && b <= r) {
      delay[k] = composite(delay[k], v);
      propagate(k, h);
      return;
    }
    propagate(k, h);
    if (b <= l || r <= a) {
      return;
    }  //
    update(l, r, v, 2 * k, a, (a + b) / 2, h - 1);
    update(l, r, v, 2 * k + 1, (a + b) / 2, b, h - 1);
    val[k] = calc(val[2 * k], val[2 * k + 1]);
  }
  inline Val query(int l, int r) { return query(l, r, 1, 0, n, height); }
  Val query(int l, int r, int k, int a, int b, int h) {
    if (b <= l || r <= a) {
      return calc.id();
    }  //
    propagate(k, h);
    if (l <= a && b <= r) {
      return val[k];
    }
    Val vall = query(l, r, 2 * k, a, (a + b) / 2, h - 1);
    Val valr = query(l, r, 2 * k + 1, (a + b) / 2, b, h - 1);
    return calc(vall, valr);
  }
};
struct Edge {
  int to, cap;
  LL cost;
  int rev;
};
class Solver {
 public:
  int n;
  vector<LL> a, x, y;
  LazySegmentTree<Min_Add> seg;  //解説中にあるSegment Tree C,Dタイプの最小コストパスを求めるために使う
  LazySegmentTree<Min_Set> segA, segB; //A,Bタイプ
  multiset<int> Qpos; //残余グラフ上の$ Q\_i \to T$容量と同じ個数の i を含む多重集合
  Solver(LL n) : n(n), a(n), x(n), y(n), seg(0), segA(0), segB(0){};

  void solve() {
    input();
    makeGraph();
    //vertex
    //P:0~n-1 (P_i = i)
    //Q:n~2*n (Q_i = n + i)
    //S:2*n+1,T:2*n+2

    //path
    //A:S->P->Q->Q->T
    //B:S->P->Q_{n}->Q->T
    //C:S->P->Q_i->Q_j->P->Q_{n}->Q->T (i>j)
    //D:S->P->Q_i->Q_j->P->Q_{n}->Q->T (i<j)

    LL mincost = 0;
    for (int i = 1; i <= n; ++i) {
      auto path = minCostAugment();
      mincost += augmentPath(path);
      cout << -mincost << endl;
    }
  }
  vector<int> minCostAugment() {
    int idA = 0, idB = 0, idCl = 0, idCr = 0, idDl = 0, idDr = 0;
    LL costA = LINF, costB = LINF, costC = LINF, costD = LINF;
    tie(costA, idA) = segA.query(*Qpos.begin(), n);
    tie(costB, idB) = segB.query(0, n);
    Node node = seg.query(0, n + 1);
    tie(costC, idCl, idCr) = node.Dc.ref();
    if (node.E.val == 0) {
      tie(costD, idDl, idDr) = node.CdNoMinE.ref();
    } else {
      tie(costD, idDl, idDr) = min(node.CdMinE, node.CdNoMinE).ref();
    }
    vector<int> ret;
    if (min({costA, costB, costC, costD}) == costA) {
      ret = {2 * n + 1, idA, n + idA};
      //DBG('A', ret)
    } else if (min({costA, costB, costC, costD}) == costB) {
      ret = {2 * n + 1, idB, n + n};
      //DBG('B', ret)
    } else if (min({costA, costB, costC, costD}) == costC) {
      ret = {2 * n + 1, idCr, idCr + n, idCl + n, idCl, n + n};
      //DBG('C', ret)
    } else {
      ret = {2 * n + 1, idDl, idDl + n, idDr + n, idDr, n + n};
      //DBG('D', ret)
    }
    connecttoT(ret);
    return ret;
  }
  LL augmentPath(vector<int> &p) {
    LL ret = 0;
    for (int i = 0; i < (int)p.size() - 1; ++i) {
      int u = p[i], v = p[i + 1];
      if (u == 2 * n + 1) {  //S->P
        segA.update(v, v + 1, LINF);
        segB.update(v, v + 1, LINF);
        seg.update(v, v + 1, Operator(x[v] + LINF, 0, 0));
      } else if (v == 2 * n + 2) {  //Q->T
        Qpos.erase(Qpos.lower_bound(u - n));
      } else if (0 <= u && u < n && v == n + n) {  //P -> Q_n
        ret -= y[u];
      } else if (0 <= u && u < n && n <= v && v < n + n) {  //P->Q
        seg.update(u, u + 1, Operator(0, -LINF / 2 + x[u] - y[u], 0));
        ret -= x[u];
      } else if (0 <= v && v < n && n <= u && u < n + n) {  //Q->P
        seg.update(v, v + 1, Operator(0, +LINF / 2 - x[v] + y[v], 0));
        ret += x[v];
      } else if (n <= u && u <= n + n && n <= v && v <= n + n && u < v) {  //Q_i->Q_j (i<j)
        seg.update(u - n, v - n, Operator(0, 0, -1));
      } else if (n <= u && u <= n + n && n <= v && v <= n + n && u > v) {  //Q_i->Q_j (i>j)
        seg.update(v - n, u - n, Operator(0, 0, 1));
      } else if (n <= u && u <= n + n && n <= v && v <= n + n && u == v) {  //Q_i->Q_i
      } else {
        assert(false);
      }
    }
    return ret;
  }

  void connecttoT(vector<int> &v) {
    auto i = *prev(Qpos.upper_bound(v.back() - n));
    v.push_back(i + n);
    v.push_back(2 * n + 2);
  }
  void input() {
    vector<vector<LL>> tmp(n, vector<LL>(3));
    for (int i = 0; i < n; ++i) {
      cin >> tmp[i][0] >> tmp[i][1] >> tmp[i][2];
    }
    sort(tmp.begin(), tmp.end());
    for (int i = 0; i < n; ++i) {
      a[i] = tmp[i][0] - 1;
      x[i] = tmp[i][1];
      y[i] = tmp[i][2];
    }
  }

  void makeGraph() {
    for (int i = 0; i < n; ++i) {
      auto lb = distance(a.begin(), lower_bound(a.begin(), a.end(), i));
      Qpos.insert(lb);
    }
    vector<Node> nodes(n);
    vector<pair<LL, int>> costA(n), costB(n);
    for (int i = 0; i < n; ++i) {
      nodes[i].E = Min_Range(0, i, i);
      nodes[i].CL = Min_Id(-x[i], i);
      nodes[i].CR = Min_Id(LINF, -1);
      nodes[i].DL = Min_Id(LINF / 2, i);
      nodes[i].DR = Min_Id(LINF, -1);
      nodes[i].Dc = Min_IdLR(LINF * 2, -1, -1);
      nodes[i].CdMinE = Min_IdLR(LINF * 2, -1, -1);
      nodes[i].CdNoMinE = Min_IdLR(LINF * 2, -1, -1);
      costA[i] = {-x[i], i};
      costB[i] = {-y[i], i};
    }
    seg = LazySegmentTree<Min_Add>(nodes);
    segA = LazySegmentTree<Min_Set>(costA);
    segB = LazySegmentTree<Min_Set>(costB);
  }
};

}
}  // namespace Problem

int main() {
  std::cin.tie(0);
  std::ios_base::sync_with_stdio(false);
  long long n = 0;
  std::cin >> n;

  Problem::Solver sol(n);
  sol.solve();
  return 0;
}

note

新たなモノイドが手に入った。