常用算法模板

大量早期代码，码风很差 ......

有些东西放这儿只是为了 CF 复制方便。

倍增后缀数组，SA Doubling

这里贴了新的板子，包含了求 height(LCP) 和 ST 表维护任意 LCP 的过程，稍微写一下主要的问题。

\(sa_i\) 是第 \(i\) 名的后缀，\(rk_i\) 是 \({\rm suffix}(i)\) 的排名。

radixSort() 是基数排序过程，第一关键字是 \(rk\)，第二关键字是 \(tp\)，其中 \(tp_i\) 表示第二关键字第 \(i\) 大的后缀是谁，排序时只需要按第一关键字做前缀和，第二关键字从后往前加入即可。

先按字符排序，随后倍增。按照当前 SA 的顺序加入第二关键字，排序后计算新的 rank。注意第一和第二关键字都相同算作一样，如果两者中任一或全部都没有第二关键字，视为不同，因为长度不同。

计算 height 就十分简单了，按照原串中的顺序依次计算，性质显然。

为了使用 ST 表，需要全局预处理出 \(\log_2 x\)。

复杂度 \(\mathcal O(n \log n)\)

int lg[N];
void preLog2(){
  for(int i=2; i<Len; i<<=1) lg[i] = 1;
  for(int i=3; i<Len; ++i) lg[i] += lg[i - 1];
}

struct SuffixArray{
  char s[Len];
  int sa[Len], rk[Len], ht[Len], st[lgLen][Len], tp[Len];
  int n, m;

  void radixSort(){
    static int pos[N];
    fill_n(pos, m + 1, 0);
    for(int i=0; i<n; i++) ++pos[rk[i]];
    for(int i=1; i<=m; i++) pos[i] += pos[i - 1];
    for(int i=n-1; i>=0; i--) sa[--pos[rk[tp[i]]]] = tp[i];
  }
  void buildSA(){
    m = 26;
    for(int i=0; i<n; i++) rk[i] = s[i] - 'a', tp[i] = i;
    radixSort();

    for(int w = 1, p = 0; p + 1 < n; m = p, w <<= 1){
      p = 0;
      for(int i=n-w; i<n; i++) tp[p++] = i;
      for(int i=0; i<n; i++)
        if(sa[i] >= w) tp[p++] = sa[i] - w;
      radixSort();

      copy_n(rk, n, tp);
      rk[sa[0]] = 0; p = 0;
      for(int i=1; i<n; i++)
        if(tp[sa[i]] == tp[sa[i - 1]] &&
        sa[i] + w < n && sa[i-1] + w < n &&
        tp[sa[i] + w] == tp[sa[i - 1] + w])
          rk[sa[i]] = p;
        else rk[sa[i]] = ++p;
    }

    int p = 0;
    for(int i=0; i<n; i++){
      if(p) --p;
      if(rk[i] == 0) continue;
      int j = sa[rk[i] - 1];
      while(s[i + p] == s[j + p]) ++p;
      ht[rk[i]] = p;
    }
  }
  void buildST(){
    for(int i=0; i<n; i++) st[0][i] = ht[i];
    for(int i=1; i<=lg[n]; i++)
      for(int j=0; j+(1<<i)-1 < n; j++)
        st[i][j] = min(st[i-1][j], st[i-1][j + (1 << (i-1))]);
  }
  int lcp(int x, int y){
    if(x == y) return n - x;
    x = rk[x], y = rk[y];
    if(x > y) swap(x, y);
    ++x;
    int k = lg[y - x + 1];
    return min(st[k][x], st[k][y - (1 << k) + 1]);
  }
}

exgcd

int exgcd(int a, int b, int &x, int &y){
    if(b == 0){
        x = 1; y = 0;
        return a;
    }
    int result = exgcd(b, a%b, x, y);
    int tp = x; x = y; y = tp - a/b*y;
    return result;
}

tarjan 缩点

void tarjan(int x){
  low[x] = dfn[x] = ++id;
  stk[++top] = x; instk[x] = true;
  for(int i=hd[x]; i; i=nt[i]){
    if(dfn[to[i]] == 0) tarjan(to[i]), low[x] = min(low[x], low[to[i]]);
    else if(instk[to[i]]) low[x] = min(low[x], dfn[to[i]]);
   }
   if(low[x] == dfn[x]){
    ++color; col[x] = color;
    while(stk[top] != x) instk[stk[top]] = false, col[stk[top--]] = color;
    top--; instk[x] = false;
  }
}

tarjan 割点/点双连通分量

以下代码统计每个 BCC 中的割点数量，摘自题目Luogu P3225 [HNOI2012]矿场搭建

void tarjan(int x, int fa){
  dfn[x] = low[x] = ++ID;
  stk[++top] = x;
  int ch = 0, last = 0;
  for(int i=hd[x]; i; i=nt[i]){
    int v = to[i];
    if(!dfn[v]){
      ++ch;
      tarjan(v, x);
      low[x] = min(low[x], low[v]);
      if(x == fa && ch >= 2){
        cut[x] = true;
        if(ch == 2 && last != 0) ++cutCnt[last]; //是根且有2个或以上孩子，是割点
      }
      if(low[v] >= dfn[x]){ //是根或没有返回边，得到一个点双
        if(x != fa) cut[x] = true; //不是根就是割点
        ++num; last = num; cnt[num] = cutCnt[num] = 0;
        while(stk[top] != v) cutCnt[num] += (int)(cut[stk[top--]]), ++cnt[num];
        cutCnt[num] += (int)(cut[stk[top--]]); ++cnt[num]; //注意加入v
        cutCnt[num] += (int)cut[x]; ++cnt[num]; //注意加入自己
      }
    }
    else if(v != fa)
      low[x] = min(low[x], dfn[v]);
  }
}

ISAP 最大流

template<class cap_t, class val_t, size_t node_c, size_t edge_c>
struct maxflow_graph {
  static_assert(sizeof(val_t) >= sizeof(cap_t), "result may overflow");

  struct edge {
    int v;
    cap_t f;
  } E[edge_c * 2];
  vector<int> G[node_c];
  size_t cur[node_c];
  int ec = 0;
  int s, t, nc;
  int h[node_c], gap[node_c + 1];

  void clear() {
    ec = 0;
    memset(cur, 0, sizeof(cur));
    for (size_t i = 0; i < node_c; ++i) G[i].clear();
  }

  void link(int x, int y, cap_t flow) {
    if (!flow) return;
    G[x].push_back(ec); E[ec++] = {y, flow};
    G[y].push_back(ec); E[ec++] = {x, 0};
  }

  bool bfs() {
    memset(h, -1, sizeof(int) * nc);
    h[t] = 0;
    int *q = new int[node_c], *ql = q - 1, *qr = q;
    *qr = t;
    while (ql != qr) {
      int x = *(++ql);
      ++gap[h[x]];
      for (int i : G[x]) {
        if (h[E[i].v] == -1 && !E[i].f)
          h[E[i].v] = h[x] + 1, *(++qr) = E[i].v;
      }
    }
    delete[] q;
    return h[s] != -1;
  }

  val_t augment(int x, val_t flow) {
    if (x == t) return flow;
    val_t res = 0;
    for (size_t &i = cur[x]; i != G[x].size(); ++i) {
      edge &e = E[G[x][i]];
      if (!e.f || h[e.v] + 1 != h[x]) continue;
      cap_t aug = min<val_t>(e.f, flow);
      aug = augment(e.v, aug);
      res += aug; flow -= aug;
      e.f -= aug; E[G[x][i] ^ 1].f += aug;
      if (!flow) return res;
    }
    --gap[h[x]];
    if (!gap[h[x]]) h[s] = nc;
    ++gap[++h[x]];
    cur[x] = 0;
    return res;
  }

  val_t maxflow() {
    val_t flow = 0;
    if (!bfs()) return 0;
    while (h[s] < nc) flow += augment(s, numeric_limits<val_t>::max());
    return flow;
  }
};

最小费用最大流

又长又慢

建议写 Dijkstra 版，不要用这个！

struct MinCostMaxFlow{
  int to[Edge], nt[Edge], hd[Node], f[Edge], w[Edge], h[Node], cur[Node], fromE[Node], ec = 0;
  int S, T;
  void clear(){memset(hd, -1, sizeof(hd)); ec = 0;}
  MinCostMaxFlow(){clear();}
  void link(int x, int y, int flow, int cost){
    to[ec] = y; nt[ec] = hd[x]; hd[x] = ec; f[ec] = flow; w[ec] = cost; ++ec;
    to[ec] = x; nt[ec] = hd[y]; hd[y] = ec; f[ec] = 0; w[ec] = -cost; ++ec;
  }

  bool spfa(){
    deque<int> q;
    memset(h, 0x3f, sizeof(h));
    h[S] = 0; q.push_back(S);
    while(!q.empty()){
      int x = q.front(); q.pop_front();
      for(int i=hd[x]; i!=-1; i=nt[i]){
        if(f[i] && h[to[i]] > h[x] + w[i]){
          h[to[i]] = h[x] + w[i];
          fromE[to[i]] = i;
          if(!q.empty() && h[to[i]] <= h[q.front()]) q.push_front(to[i]);
          else q.push_back(to[i]);
        }
      }
    }
    return h[T] != 0x3f3f3f3f;
  }
  bool vis[Node];
  int dfs(int x, int flow){
    if(x == T || !flow) return flow;
    int res = 0;
    vis[x] = true;
    for(int i=hd[x]; i!=-1; i=nt[i]){
      if(vis[to[i]] || h[to[i]] != h[x] + w[i] || !f[i]) continue;
      int aug = min(flow, f[i]);
      if((aug = min(aug, dfs(to[i], aug)))){
        f[i] -= aug; f[i^1] += aug;
        flow -= aug; res += aug;
      }
      if(!flow) break;
    }
    vis[x] = false;
    if(!res) h[x] = 0x3f3f3f3f;
    return res;
  }
  pair<int, int> mcmf(){
    int flow = 0, cost = 0;
    while(spfa()){
      int fl = dfs(S, 0x3f3f3f3f);
      flow += fl; cost += fl * h[T];
    }
    return make_pair(flow, cost);
  }
}G;

回文自动机 PAM

namespace pam{
char s[N];
int c[N][26], fa[N], len[N], cnt[N], last, pos = 0, nc = 0;

void init() {
  len[0] = 0, len[1] = -1;
  fa[0] = 1;
  s[0] = -1;
  nc = 1;
  last = 0;
}
int getfail(int x) {
  while (s[pos] != s[pos - len[x] - 1])
    x = fa[x];
  return x;
}

void extend(int x) {
  s[++pos] = x;
  int now = getfail(last);
  if(!c[now][x]) {
    ++nc;
    len[nc] = len[now] + 2;
    fa[nc] = c[getfail(fa[now])][x];
    c[now][x] = nc;
  }
  last = c[now][x];
  ++cnt[last];
}
}; // namespace pam

模素数组合数 / 二项式系数

int fac[N], ifac[N];
void prefac(int n) {
  fac[0] = 1;
  for (int i = 1; i <= n; ++i) fac[i] = mul(fac[i - 1], i);
  ifac[n] = inv(fac[n]);
  for (int i = n - 1; i >= 0; --i) ifac[i] = mul(ifac[i + 1], i + 1);
}
int binom(int n, int m) {
  return mul(fac[n], ifac[m], ifac[n - m]);
}

SAIS 线性后缀数组

封装了求 LCP，但注意 sais() 函数需要保证字符串是 int 类型，且最后存在极小字符，最后极小字符被排在开头。

（即 0 下标变 1 下标）

struct suffix_array {
  int sa[N], pos[N];
  int st[logN][N];
  const int *ht = st[0];

  void sais(const int *a, int *sa, int n, int w) {
    bool *type = new bool[n];
    int *buc = new int[w], *pl = new int[w], *pr = new int[w];
    int n1 = 0, *p1;
    type[n - 1] = true;
    for (int i = n - 1; i; --i) {
      type[i - 1] = a[i - 1] == a[i] ? type[i] : a[i - 1] < a[i];
      n1 += type[i] && !type[i - 1];
    }
    p1 = new int[n1];
    for (int i = 1, *p = p1; i < n; ++i)
      if (type[i] && !type[i - 1]) *(p++) = i;

    memset(buc, 0, sizeof(int) * w);
    for (int i = 0; i < n; ++i)
      ++buc[a[i]];
    for (int i = 0, sum = 0; i < w; ++i) {
      sum += buc[i];
      pr[i] = sum; pl[i] = sum - buc[i];
    }

    auto induce = [&]() {
      memcpy(buc, pl, sizeof(int) * w);
      for (int i = 0; i < n; ++i)
        if (int j = sa[i] - 1; j >= 0 && !type[j]) sa[buc[a[j]]++] = j;
      memcpy(buc, pr, sizeof(int) * w);
      for (int i = n - 1; i >= 0; --i)
        if (int j = sa[i] - 1; j >= 0 && type[j]) sa[--buc[a[j]]] = j;
    };

    memset(sa, -1, sizeof(int) * n);
    memcpy(buc, pr, sizeof(int) * w);
    for (int i = 1; i < n; ++i)
      if (type[i] && !type[i - 1]) sa[--buc[a[i]]] = i;
    induce();

    int id = 0;
    int *s1 = new int[n1], *sa1 = new int[n1];
    for (int i = 0, *p = sa1; i < n; ++i)
      if (type[sa[i]] && sa[i] && !type[sa[i] - 1]) *(p++) = sa[i];
    sa[sa1[0]] = 0;
    for (int i = 0; i < n1 - 1; ++i) {
      int p = sa1[i], q = sa1[i + 1];
      if (a[p] != a[q]) {
        sa[q] = ++id;
        continue;
      }
      ++p, ++q;
      while (a[p] == a[q] && type[p] == type[q] && !(type[p] && !type[p - 1]) && !(type[q] && !type[q - 1]))
        ++p, ++q;
      id += !((type[p] && !type[p - 1]) && (type[q] && !type[q - 1]));
      sa[sa1[i + 1]] = id;
    }
    for (int i = 0; i < n1; ++i)
      s1[i] = sa[p1[i]];

    if (id + 1 == n1)
      for (int i = 0; i < n1; ++i)
        sa1[s1[i]] = i;
    else
      sais(s1, sa1, n1, id + 1);

    memcpy(buc, pr, sizeof(int) * w);
    memset(sa, -1, sizeof(int) * n);
    for (int i = n1 - 1; i >= 0; --i)
      sa[--buc[a[p1[sa1[i]]]]] = p1[sa1[i]];
    induce();
    delete[] type; delete[] buc; delete[] pl; delete[] pr; delete[] p1; delete[] s1; delete[] sa1;
  }

  void init(char *s) {
    int n = strlen(s);
    int *a = new int[n + 1];
    copy(s, s + n, a);
    a[n] = 0;
    sais(a, sa, n + 1, 128);
    delete[] a;

    for (int i = 1; i <= n; ++i)
      pos[sa[i]] = i;
    int h = 0;
    for (int i = 0; i < n; ++i) {
      h = !h ? 0 : h - 1;
      int j = sa[pos[i] - 1];
      while (s[i + h] == s[j + h]) ++h;
      st[0][pos[i]] = h;
    }

    for (int i = 1, li = __lg(n); i <= li; ++i)
      for (int j = 1, li = n - (1 << i) + 1; j <= li; ++j)
        st[i][j] = min(st[i - 1][j], st[i - 1][j + (1 << (i - 1))]);
  }

  int lcp(int i, int j) {
    i = pos[i], j = pos[j];
    if (i > j) swap(i, j);
    int k = __lg(j - i);
    return min(st[k][i + 1], st[k][j - (1 << k) + 1]);
  }
}

复数 FFT

using cp = complex<double>;
const double pi = acos(-1);

cp root[N];

cp MUL(const cp &a, const cp &b) {
  double x = a.real(), y = a.imag(), u = b.real(), v = b.imag();
  return cp(x * u - y * v, x * v + y * u);
}

void init_roots() {
  for (int i = N / 2; i; i >>= 1) {
    double theta = pi / i;
    for (int j = 0; j < i; ++j)
      root[i + j] = polar(1.0, j * theta);
  }
}
void DIF(cp *a, int n) {
  for (int l = n, k = l / 2; k; l >>= 1, k >>= 1) {
    cp *w = root + k;
    for (cp *i = a, *t = a + n; i != t; i += l)
      for (int j = 0; j < k; ++j) {
        cp t = i[j];
        i[j] += i[j + k];
        i[j + k] = MUL((t - i[j + k]), w[j]);
      }
  }
}
void DIT(cp *a, int n) {
  for (int l = 2, k = 1; k < n; l <<= 1, k <<= 1) {
    cp *w = root + k;
    for (cp *i = a, *t = a + n; i != t; i += l)
      for (int j = 0; j < k; ++j) {
        cp t = MUL(i[j + k], w[j]);
        i[j + k] = i[j] - t;
        i[j] += t;
      }
  }
  double iv = 1.0 / n;
  reverse(a + 1, a + n);
  for (int i = 0; i < n; ++i)
    a[i] *= iv;
}

待补充...