# C++S's CODEBOOK
###### by @XDEv11, @D4nnyLee
[TOC]
```
"vimrc
if filereadable("/etc/vim/vimrc")
source /etc/vim/vimrc
endif
colo default
set bg=dark
set cul
hi cursorline cterm=bold ctermbg=none ctermfg=none
set cuc
hi cursorcolumn cterm=bold ctermbg=none ctermfg=none
set is
set hls
hi search cterm=underline ctermbg=none ctermfg=none
set nu
set cindent
set aw
set ts=4
set shiftwidth=4
set sts=0
set noet
imap <F1> <C-o>ofor (int i{0}; i < ; ++i)<C-o>F;
imap <F2> <C-o>ofor (auto& x : )<left>
imap jk <ESC>
```
```
## makefile
CXXFLAGS := -std=c++17 -O2 -Wall -Wextra -Wshadow \
--debug -fsanitize=address,undefined -D_GLIBCXX_DEBUG
MAKE_FLAG := --no-print-directory
%: INPUT = $(shell echo $@ | tr [:upper:] [:lower:]).in
%: EXE = x$@.exe
%:
@make $(MAKE_FLAG) $(EXE)
@make $(MAKE_FLAG) $(INPUT)
./$(EXE) < $(INPUT)
x%.exe: %.cpp
g++ $(CXXFLAGS) -o $@ $<
%.in:
vim $@
clean:
rm -f *.exe
```
```bash!
#!/usr/bin/env bash
CODE=$1.cpp
INPUT=`echo $1 | tr [:upper:] [:lower:]`.in
EXE=x$1.exe
CXXFLAGS='-std=c++17 -O2 -Wall -Wextra -Wshadow
--debug -fsanitize=address,undefined -D_GLIBCXX_DEBUG'
g++ $CXXFLAGS -o $EXE $CODE
if [[ $? -eq 0 ]]; then
if [[ ! -e $INPUT ]]; then
vim $INPUT
fi
echo -e 'Executing...\n'
./$EXE < $INPUT
fi
```
`#include <bits/stdc++.h>`
# Data Structure
### DSU
```cpp!
// fast disjoint set union
class DSU {
vector<int> pa, sz;
public:
explicit DSU(int n) : pa(n, -1), sz(n, 1) {}
int find(int x) { // collapsing find
return pa[x] == -1 ? x : pa[x] = find(pa[x]);
}
void unionn(int x, int y) { // weighted union
auto rx{find(x)}, ry{find(y)};
if (rx == ry) return ;
if (sz[rx] < sz[ry]) swap(rx, ry);
pa[ry] = rx, sz[rx] += sz[ry];
}
};
```
```cpp!
// previous/next one
class PvNx {
vector<int> pa, sz, mn, mx;
int find(int x) { // collapsing find
return pa[x] == -1 ? x : pa[x] = find(pa[x]);
}
void unionn(int x, int y) { // weighted union
auto rx{find(x)}, ry{find(y)};
if (rx == ry) return ;
if (sz[rx] < sz[ry]) swap(rx, ry);
pa[ry] = rx, sz[rx] += sz[ry], mn[rx] = min(mn[rx], mn[ry]), mx[rx] = max(mx[rx], mx[ry]);
}
public:
explicit PvNx(int n) : pa(n + 1, -1), sz(n + 1, 1), mn(n + 1) { iota(mn.begin(), mn.end(), 0), mx = mn; }
void remove(int i) { unionn(i, i + 1); }
int prev(int i) { return mn[find(i)] - 1; }
int next(int i) {
int j{mx[find(i)]};
if (i == j) j = mx[find(j + 1)];
return j;
}
bool exist(int i) { return i == mx[find(i)]; }
};
```
### Sparse Table
```cpp!
// sparse table
template<typename value_t, class merge_t>
class ST {
int n;
vector<int> log2;
vector<vector<value_t>> t;
merge_t merge; // associative & idempotent function for ST
public:
explicit ST(const vector<value_t>& v, const merge_t& _merge = merge_t{})
: n{v.size()}, log2(n + 1), merge{_merge} {
log2[1] = 0;
for (int i{2}; i <= n; ++i) log2[i] = log2[i / 2] + 1;
t.assign(n, vector<value_t>(log2[n] + 1));
for (int i{n - 1}; i >= 0; --i) {
t[i][0] = v[i];
for (int j{1}; i + (1 << j) <= n; ++j)
t[i][j] = merge(t[i][j - 1], t[i + (1 << (j - 1))][j - 1]);
}
}
value_t query(int l, int r) { // [l:r)
int j{log2[r - l]};
return merge(t[l][j], t[r - (1 << j)][j]);
}
};
```
### Segment Tree
> [source](https://codeforces.com/blog/entry/18051)
> [proof](https://hackmd.io/@XDEv11/iterative-segment-tree-proof)
```cpp!
// segment tree
// segment tree
template<typename value_t, class merge_t>
class SGT {
int n;
vector<value_t> t; // root starts at 1
value_t defa;
merge_t merge; // associative function for SGT
public:
explicit SGT(int _n, value_t _defa, const merge_t& _merge = merge_t{})
: n{_n}, t(2 * n), defa{_defa}, merge{_merge} {}
void modify(int p, const value_t& x) {
for (t[p += n] = x; p > 1; p >>= 1) t[p >> 1] = merge(t[p], t[p ^ 1]);
}
value_t query(int l, int r) { return query(l, r, defa); }
value_t query(int l, int r, value_t init) { // [l:r)
for (l += n, r += n; l < r; l >>= 1, r >>= 1) {
if (l & 1) init = merge(init, t[l++]);
if (r & 1) init = merge(init, t[--r]);
}
return init;
}
};
```
```cpp!
class SGT {
int n;
vector<long long> t;
int left(int tv) { return tv + 1; }
int right(int tv, int tl, int tm) { return tv + 2 * (tm - tl); }
void modify(int p, long long x, int tv, int tl, int tr) {
if (tr - tl > 1) {
int tm{(tl + tr) / 2};
if (p < tm) modify(p, x, left(tv), tl, tm);
else modify(p, x, right(tv, tl, tm), tm, tr);
t[tv] = t[left(tv)] + t[right(tv, tl, tm)];
} else t[tv] = x;
}
long long query(int l, int r, int tv, int tl, int tr) {
if (l == tl && r == tr) return t[tv];
int tm{(tl + tr) / 2};
if (r <= tm) return query(l, r, left(tv), tl, tm);
else if (l >= tm) return query(l, r, right(tv, tl, tm), tm, tr);
else return query(l, tm, left(tv), tl, tm) +
query(tm, r, right(tv, tl, tm), tm, tr);
}
public:
explicit SGT(int _n) : n{_n}, t(2 * n - 1) {}
void modify(int p, long long x) { modify(p, x, 0, 0, n); };
long long query(int l, int r) { return query(l, r, 0, 0, n); }
int ps_lower_bound(long long ps) { // prefix sum lower bound
if (ps > t[0]) return n;
int tv{0}, tl{0}, tr{n};
while (tr - tl > 1) {
int tm{(tl + tr) / 2};
if (t[left(tv)] >= ps) tv = left(tv), tr = tm;
else ps -= t[left(tv)], tv = right(tv, tl, tm), tl = tm;
}
return tl;
}
};
```
#### range modification (lazy propagation)
```cpp!
// segment tree
// range query & range modify
template<typename value_t>
class SGT {
using node_t = pair<value_t, int>;
int n;
vector<node_t> t;
vector<optional<value_t>> lz;
// [ tv+1 : tv+2*(tm-tl) ) -> left subtree
int left(int tv) { return tv + 1; }
int right(int tv, int tl, int tm) { return tv + 2 * (tm - tl); }
/** differ from case to case **/
// query is "max" and modify is "add" here
node_t merge(const node_t& x, const node_t& y) { // associative function
return max(x, y);
}
void update(int tv, int len, const value_t& x) {
if (!lz[tv]) lz[tv] = x;
else lz[tv] = lz[tv].value() + x;
t[tv].fi = t[tv].fi + x;
}
/******************************/
void build(const vector<value_t>& v, int tv, int tl, int tr) {
if (tr - tl > 1) {
int tm{(tl + tr) / 2};
build(v, left(tv), tl, tm);
build(v, right(tv, tl, tm), tm, tr);
t[tv] = merge(t[left(tv)], t[right(tv, tl, tm)]);
} else t[tv] = {v[tl], tl};
}
void push(int tv, int tl, int tr) { // lazy propagation
if (!lz[tv]) return ;
int tm{(tl + tr) / 2};
update(left(tv), tm - tl, lz[tv].value());
update(right(tv, tl, tm), tr - tm, lz[tv].value());
lz[tv].reset();
}
void set(int p, const value_t& x, int tv, int tl, int tr) {
if (tr - tl > 1) {
push(tv, tl, tr);
int tm{(tl + tr) / 2};
if (p < tm) set(p, x, left(tv), tl, tm);
else set(p, x, right(tv, tl, tm), tm, tr);
t[tv] = merge(t[left(tv)], t[right(tv, tl, tm)]);
} else t[tv].fi = x;
}
void rmodify(int l, int r, const value_t& x, int tv, int tl, int tr) {
if (!(l == tl && r == tr)) {
push(tv, tl, tr);
int tm{(tl + tr) / 2};
if (r <= tm) rmodify(l, r, x, left(tv), tl, tm);
else if (l >= tm) rmodify(l, r, x, right(tv, tl, tm), tm, tr);
else rmodify(l, tm, x, left(tv), tl, tm),
rmodify(tm, r, x, right(tv, tl, tm), tm, tr);
t[tv] = merge(t[left(tv)], t[right(tv, tl, tm)]);
} else update(tv, tr - tl, x);
}
node_t rquery(int l, int r, int tv, int tl, int tr) {
if (l == tl && r == tr) return t[tv];
push(tv, tl, tr);
int tm{(tl + tr) / 2};
if (r <= tm) return rquery(l, r, left(tv), tl, tm);
else if (l >= tm) return rquery(l, r, right(tv, tl, tm), tm, tr);
else return merge(rquery(l, tm, left(tv), tl, tm),
rquery(tm, r, right(tv, tl, tm), tm, tr));
}
public:
explicit SGT(const vector<value_t>& v) : n{v.size()}, t(2 * n - 1), lz(2 * n - 1) { build(v, 0, 0, n); }
void set(int p, const value_t& x) { set(p, x, 0, 0, n); }
void rmodify(int l, int r, const value_t& x) { rmodify(l, r, x, 0, 0, n); } // [l:r)
node_t rquery(int l, int r) { return rquery(l, r, 0, 0, n); } // [l:r)
};
```
### Time Segment Tree (with Intermittent Event)
```cpp!
class DSU {
vector<int> pa, sz;
stack<pair<int, int>> rec;
public:
explicit DSU(int n) : pa(n, -1), sz(n, 1), rec{} {}
int find(int x) {
return pa[x] == -1 ? x : find(pa[x]);
}
bool unionn(int x, int y) { // weighted union
auto rx{find(x)}, ry{find(y)};
if (rx == ry) return false;
if (sz[rx] < sz[ry]) swap(rx, ry);
pa[ry] = rx, sz[rx] += sz[ry];
rec.emplace(rx, ry);
return true;
}
void undo() {
auto [rx, ry]{rec.top()}; rec.pop();
pa[ry] = -1, sz[rx] -= sz[ry];
}
};
class EventSGT {
int n, m;
vector<vector<pair<int, int>>> t;
DSU dsu;
int left(int tv) { return tv + 1; }
int right(int tv, int tl, int tm) { return tv + 2 * (tm - tl); }
void add_event(pair<int, int> e, int l, int r, int tv, int tl, int tr) {
if (l == tl && r == tr) return t[tv].push_back(e), []{}();
int tm{(tl + tr) / 2};
if (r <= tm) add_event(e, l, r, left(tv), tl, tm);
else if (l >= tm) add_event(e, l, r, right(tv, tl, tm), tm, tr);
else add_event(e, l, tm, left(tv), tl, tm), add_event(e, tm, r, right(tv, tl, tm), tm, tr);
}
void run(int tv, int tl, int tr) {
int cnt{};
for (auto& [x, y] : t[tv]) cnt += dsu.unionn(x, y);
if (tr - tl == 1) {
auto& [x, y]{query[tl]};
if (x != -1 && y != -1) ans[tl] = (dsu.find(x) == dsu.find(y));
} else {
int tm{(tl + tr) / 2};
run(left(tv), tl, tm), run(right(tv, tl, tm), tm, tr);
}
while (cnt--) dsu.undo();
}
public:
vector<pair<int, int>> query;
vector<bool> ans;
explicit EventSGT(int _n, int _m) : n{_n}, m{_m}, t(2 * n - 1), dsu{m}, query(n, {-1, -1}), ans(n) {}
void add_event(pair<int, int> e, int l, int r) { add_event(e, l, r, 0, 0, n); }
void add_query(pair<int, int> e, int tm) { query[tm] = e; }
void run() { run(0, 0, n); }
bool get_ans(int tm) { return ans[tm]; }
};
```
# Algorithm
## Sequence
### LIS
```cpp!
// longest increasing subsequence
int LIS(const vector<int>& s) {
// use binary search to update the smallest element
// that ends in a subsequence of length d
vector<int> v{};
v.push_back(s[0]);
for (int i{1}; i < s.size(); ++i)
if (s[i] > v.back()) v.push_back(s[i]);
else *lower_bound(v.begin(), v.end(), s[i]) = s[i];
return v.size();
}
```
### LCS
```cpp!
// longest common subsequence
int LCS(const vector<int>& a, const vector<int>& b) {
vector<vector<int>> dp(a.size() + 1, vector<int>(b.size() + 1));
for (int i{1}; i <= int(a.size()); ++i)
for (int j{1}; j <= int(b.size()); ++j) {
dp[i][j] = max(dp[i][j], max(dp[i - 1][j], dp[i][j - 1]));
if (a[i - 1] == b[j - 1]) dp[i][j] = max(dp[i][j], dp[i - 1][j - 1] + 1);
}
return dp[a.size()][b.size()];
}
```
## Graph
### tree diameter
```cpp!
// tree diameter
int dfs(int u, int w = -1) {
int mx{0};
for (auto& v : adj[u])
if (v != w) {
auto len{1 + dfs(v, u)};
diam = max(diam, mx + len);
mx = max(mx, len);
}
return mx;
}
```
```cpp!
pair<int, int> dfs(int u, int w = -1) {
pair<int, int> mx{0, u};
for (auto& v : adj[u])
if (v != w) {
auto [len, leaf]{dfs(v, u)};
mx = max(mx, {len + 1, leaf});
}
return mx;
}
int tree_diameter() {
// assume 'b' and 'c' are endpoints of a tree diameter
// choose a arbitrary node 'a'
// let 'x' be the node that 'a' first meet path "b - c"
// 'a' will finally find either 'b' or 'c' through 'x' as its longest path
auto b{dfs(0).second};
return dfs(b).first;
}
```
### all longest path
```cpp!
// all longest path (generalization of the tree diameter problem)
vector<tuple<int, int, int>> dp{};
// [mx1, x, mx2] the path of mx1 goes through x
int dfs1(int u, int w = -1) {
int mx{0};
for (auto& v : adj[u])
if (v != w) {
auto l{1 + dfs1(v, u)};
mx = max(mx, l);
auto& [mx1, x, mx2]{dp[u]};
if (l >= mx1) mx2 = mx1, mx1 = l, x = v;
else if (l > mx2) mx2 = l;
}
return mx;
}
void dfs2(int u, int w = -1) {
if (w != -1) {
int tmx;
{ // calculate the longest path through parent
auto& [mx1, x, mx2]{dp[w]};
if (x != u) tmx = mx1 + 1;
else tmx = mx2 + 1;
}
{ // update the path
auto& [mx1, x, mx2]{dp[u]};
if (tmx >= mx1) mx2 = mx1, mx1 = tmx, x = w;
else if (tmx > mx2) mx2 = tmx;
}
}
for (auto& v : adj[u])
if (v != w) dfs2(v, u);
}
void all_longest_path() {
dfs1(0), dfs2(0);
}
```
### LCA (TODO)
```cpp!
// binary lifting
class LCA {
const vector<vector<int>>& adj;
int n;
vector<int> d;
vector<int> log2;
vector<vector<int>> an{};
void dfs(int u, int w = -1, int dep = 0) {
d[u] = dep;
an[u][0] = w;
for (int i{1}; i <= log2[n - 1] && an[u][i - 1] != -1; ++i)
an[u][i] = an[an[u][i - 1]][i - 1]; // 走 2^(i-1) 再走 2^(i-1) 步
// 因為計算 an 會用到祖先的資訊，所以先計算再繼續往下
for (auto& v : adj[u]) {
if (v == w) continue; // parent
dfs(v, u, dep + 1);
}
}
public:
LCA(const vector<vector<int>>& _adj, int root)
: adj{_adj}, n{adj.size()}, d(n), log2(n) {
log2[1] = 0;
for (int i{2}; i < log2.size(); ++i) log2[i] = log2[i / 2] + 1;
an.assign(n, vector<int>(log2[n - 1] + 1, -1));
dfs(root);
}
int operator()(int u, int v) {
if (d[u] > d[v]) swap(u, v);
for (int i{log2[d[v] - d[u]]}; i >= 0; --i)
if (d[v] - d[u] >= (1 << i)) v = an[v][i];
// v 先走到跟 u 同高度
if (u == v) return u;
for (int i{log2[d[u]]}; i >= 0; --i)
if (an[u][i] != an[v][i]) u = an[u][i], v = an[v][i];
// u, v 一起走到 lca(u, v) 的下方
return an[u][0];
// 回傳 lca(u, v)
}
};
```
```cpp!
// Euler Tour Technique
class LCA {
const vector<vector<int>>& adj;
int n;
vector<int> d, first, euler{}, log2{};
vector<vector<int>> st{};
void dfs(int u, int w = -1, int dep = 0) {
d[u] = dep;
first[u] = euler.size();
euler.push_back(u);
for (auto& v : adj[u]) {
if (v == w) continue;
dfs(v, u, dep + 1);
euler.push_back(u);
}
}
public:
LCA(const vector<vector<int>>& _adj, int root) : adj{_adj}, n{adj.size()}, d(n), first(n) {
dfs(root);
int tn{euler.size()};
log2.resize(tn + 1);
log2[1] = 0;
for (int i{2}; i <= tn; ++i) log2[i] = log2[i / 2] + 1;
st.assign(tn, vector<int>(log2[tn] + 1));
for (int i{tn - 1}; i >= 0; --i) {
st[i][0] = euler[i];
for (int j{1}; i + (1 << j) <= tn; ++j) {
auto& x{st[i][j - 1]};
auto& y{st[i + (1 << (j - 1))][j - 1]};
st[i][j] = d[x] <= d[y] ? x : y;
}
}
}
int operator()(int u, int v) {
int l{first[u]}, r{first[v]};
if (l > r) swap(l, r);
++r; // make the interval left closed right open
int j{log2[r - l]};
auto& x{st[l][j]};
auto& y{st[r - (1 << j)][j]};
return d[x] <= d[y] ? x : y;
}
};
```
```cpp!
// heavy-light decomposition
class LCA {
const vector<vector<int>>& adj;
int n;
vector<int> pa, dep, heavy, head;
// 父節點, 深度, 重邊, 路徑起點
int heavy_light(int u, int w = -1, int d = 0) {
pa[u] = w, dep[u] = d;
int mx{0}, tot{1};
for (auto& v : adj[u]) {
if (v == w) continue;
int sz{heavy_light(v, u, d + 1)};
if (sz > mx) mx = sz, heavy[u] = v;
tot += sz;
}
return tot;
}
void decompose(int u, int h) {
head[u] = h;
for (auto& v : adj[u]) {
if (v == pa[u]) continue;
if (v == heavy[u]) decompose(v, h);
else decompose(v, v);
}
}
public:
LCA(const vector<vector<int>>& _adj, int root) : adj{_adj}, n{adj.size()}, pa(n), dep(n), heavy(n), head(n) {
// Heavy-Light Decomposition
heavy_light(root);
decompose(root, root);
}
int operator()(int u, int v) {
while (head[u] != head[v]) {
if (dep[head[u]] <= dep[head[v]]) v = pa[head[v]];
else u = pa[head[u]];
}
return dep[u] <= dep[v] ? u : v;
}
};
```
### Heavy-light decomposition
```cpp!
class HLD { // Heavy-Light Decomposition
vector<int> pa, sz, dep, hv, hd;
void heavy_light(const vector<vector<int>>& adj, int u, int w = -1, int d = 0) {
pa[u] = w, sz[u] = 1, dep[u] = d;
int mx{0};
for (auto& v : adj[u]) {
if (v == w) continue;
heavy_light(adj, v, u, d + 1);
sz[u] += sz[v];
if (sz[v] > mx) mx = sz[v], hv[u] = v;
}
}
void decompose(const vector<vector<int>>& adj, int u, int h) {
hd[u] = h;
for (auto& v : adj[u]) {
if (v == pa[u]) continue;
if (v == hv[u]) decompose(adj, v, h);
else decompose(adj, v, v);
}
}
public:
HLD(const vector<vector<int>>& adj, int root) : pa(adj.size()), sz(adj.size()), dep(adj.size()), hv(adj.size(), -1), hd(adj.size()) {
heavy_light(adj, root), decompose(adj, root, root);
}
int parent(int v) { return pa[v]; }
int size(int v) { return sz[v]; }
int depth(int v) { return dep[v]; }
int heavy(int v) { return hv[v]; }
int head(int v) { return hd[v]; }
};
```
### MST
```cpp!
// Kruskal’s algorithm
vector<tuple<int, int, long long>> Kruskal(int n, vector<tuple<long long, int, int>>& edges) {
vector<tuple<int, int, long long>> mst{};
DSU dsu{n};
sort(edges.begin(), edges.end());
for (auto& [w, u, v] : edges)
if (dsu.find(u) != dsu.find(v))
dsu.unionn(u, v), mst.emplace_back(u, v, w);
return mst;
}
```
```cpp!
// Prim's algorithm
vector<tuple<int, int, long long>> Prim(const vector<vector<pair<int, long long>>>& adj) {
const auto& n{adj.size()};
vector<tuple<int, int, long long>> mst{};
vector<bool> found(n, false);
using ti = tuple<long long, int, int>;
priority_queue<ti, vector<ti>, greater<ti>> pq{};
found[0] = true;
for (auto& [v, w] : adj[0]) pq.emplace(w, 0, v);
for (int i{0}; i < n - 1; ++i) {
int mn, u, v;
do {
tie(mn, u, v) = pq.top(), pq.pop();
} while (found[v]);
found[v] = true, mst.emplace_back(u, v, mn);
for (auto& [x, w] : adj[v]) pq.emplace(w, v, x);
}
return mst;
}
```
### Articulation point & Bridge & BCC
```cpp!
// Hopcroft & Tarjan
class graph {
vector<int> dep, low;
stack<pair<int, int>> s{};
vector<int> bccno;
void dfs(int u, int w = -1, int d = 1) {
int child{0};
low[u] = dep[u] = d;
for (auto& v : adj[u]) {
if (v == w || dep[v] >= dep[u]) continue;
if (!dep[v]) { // tree edge
++child;
s.emplace(u, v);
dfs(v, u, d + 1);
low[u] = min(low[u], low[v]);
if (low[v] >= dep[u]) {
ap[u] = true;
bcc.emplace_back();
while (true) {
auto [x, y]{s.top()}; s.pop();
if (bccno[x] != bcc.size())
bcc.back().push_back(x), bccno[x] = bcc.size();
if (bccno[y] != bcc.size())
bcc.back().push_back(y), bccno[y] = bcc.size();
if (x == u && y == v) break;
}
}
if (low[v] > dep[u]) bridge.emplace_back(u, v);
} else low[u] = min(low[u], dep[v]); // back edge
}
if (w == -1) ap[u] = (child > 1);
}
public:
vector<vector<int>> adj;
vector<bool> ap;
vector<pair<int, int>> bridge{};
vector<vector<int>> bcc{};
graph(int n) : dep(n, 0), low(n), bccno(n), adj(n), ap(n, false) {}
void run() {
for (int i{0}; i < adj.size(); ++i) if (!dep[i]) dfs(i);
}
};
```
### SCC
```cpp!
// Tarjan’s strongly connected components algorithm
class SCC {
int d;
vector<int> dfn{}, low{};
vector<int> sccno{};
vector<vector<int>> scc{};
stack<int> s{};
void dfs(const vector<vector<int>>& adj, int u) {
low[u] = dfn[u] = ++d;
s.emplace(u);
for (auto& v : adj[u]) {
if (sccno[v]) continue; // cross edge
if (!dfn[v]) { // tree edge
dfs(adj, v);
low[u] = min(low[u], low[v]);
} else low[u] = min(low[u], dfn[v]); // forward / back / cross edge
}
if (low[u] == dfn[u]) { // first vertex of a scc
scc.emplace_back();
while (true) {
int x{s.top()}; s.pop();
sccno[x] = scc.size(), scc.back().push_back(x);
if (x == u) break;
}
}
}
public:
SCC()=default;
vector<vector<int>> operator()(const vector<vector<int>>& adj) {
int n{adj.size()};
d = 0, dfn.assign(n, 0), low.resize(n), sccno.assign(n, 0);
for (int i{0}; i < n; ++i) if (!dfn[i]) dfs(adj, i);
return move(scc);
}
};
```
:::spoiler Kosaraju's algorithm
```cpp!
// scc
// Kosaraju's algorithm
class graph {
vector<vector<int>> rvs; // reverse graph
vector<bool> vis;
stack<int> s{};
vector<int> component{};
void dfs1(int u) {
vis[u] = true;
for (auto& v : rvs[u])
if (!vis[v]) dfs1(v);
s.push(u);
}
void dfs2(int u) {
vis[u] = true;
component.push_back(u);
for (auto& v : adj[u])
if (!vis[v]) dfs2(v);
}
public:
vector<vector<int>> adj;
vector<int> sccno{};
vector<vector<int>> scc{};
explicit graph(int n) : rvs(n), vis(n), adj(n), sccno(n) {}
void add_edge(int u, int v) {
adj[u].push_back(v);
rvs[v].push_back(u);
}
void kosaraju() {
// find the order on reverse graph
for (int u{0}; u < adj.size(); ++u)
if (!vis[u]) dfs1(u);
// find scc on reverse order
fill(vis.begin(), vis.end(), false);
while (!s.empty()) {
auto u{s.top()}; s.pop();
if (!vis[u]) {
dfs2(u);
scc.push_back(move(component));
for (auto& v : scc.back()) sccno[v] = scc.size();
}
}
}
};
```
:::
### Topological sort
```cpp!
// topological sort 1
optional<vector<int>> top_sort(vector<vector<int>>& adj) {
vector<int> res{};
int n{static_cast<int>(adj.size())};
vector<int> cnt(n, 0); // predecessor count
for (int u{0}; u < n; ++u)
for (auto& v : adj[u]) ++cnt[v];
queue<int> qu{};
for (int u{0}; u < n; ++u) if (!cnt[u]) qu.push(u);
while (!qu.empty()) {
auto u{qu.front()}; qu.pop();
res.push_back(u);
for (auto& v : adj[u])
if (!--cnt[v]) qu.push(v);
}
if (res.size() != adj.size()) return nullopt;
return res;
}
```
```cpp!
// topological sort 2
class DAG {
public:
int n;
vector<vector<int>> adj;
vector<int> s, order{};
explicit DAG(int n) : n{n}, adj(n), s(n) {}
bool dfs(int u) {
// state 0 : unvisited, state 1 : processing, state 2 : finished
s[u] = 1;
for (auto& v : adj[u])
if (s[v] == 1 || (s[v] == 0 && !dfs(v))) return false; // cycle
order.push_back(u);
s[u] = 2;
return true;
}
bool topSort() {
for (int i{0}; i < n; ++i)
if (s[i] == 0 && !dfs(i)) return order.clear(), false;
reverse(order.begin(), order.end());
return true;
}
};
```
### 2-SAT
```cpp!
auto scc{SCC{}(adj)};
vector<int> sccid(2 * n);
for (int i{0}; i < scc.size(); ++i) {
for (auto& v : scc[i]) sccid[v] = i;
}
for (int i{0}; i < n; ++i) {
if (sccid[2 * i] == sccid[2 * i + 1]) return cout << "-1\n", true;
}
vector<vector<int>> sccadj(scc.size());
for (int i{0}; i < scc.size(); ++i)
for (auto& u : scc[i])
for (auto& v : adj[u]) {
if (sccid[v] == i) continue;
sccadj[i].push_back(sccid[v]);
}
for (auto& v : sccadj) {
sort(v.begin(), v.end());
v.erase(unique(v.begin(), v.end()), v.end());
}
vector<int> ans(n, -1);
auto ts{top_sort(sccadj)};
for (auto& i : ts.value())
for (auto& v : scc[i]) {
if (ans[v / 2] != -1) continue;
ans[v / 2] = v & 1;
}
```
:::spoiler old
```cpp!
class twoSat {
int n;
vector<vector<int>> adj;
stack<int> s{};
bool dfs(int x) {
if (mark[x ^ 1]) return false;
if (mark[x]) return true;
mark[x] = true;
s.push(x);
for (auto& y : adj[x]) if (!dfs(y)) return false;
return true;
}
public:
vector<bool> mark;
twoSat(int _n) : n{_n}, adj(2 * n), mark(2 * n, false) {}
void or_clause(int x, bool xv, int y, bool yv) { // ((x == xv) || (y == yv))
x = 2 * x + xv, y = 2 * y + yv;
adj[x ^ 1].push_back(y);
adj[y ^ 1].push_back(x);
}
void xor_clause(int x, bool xv, int y, bool yv) { // ((x == xv) xor (y == yv))
x = 2 * x + xv, y = 2 * y + yv;
adj[x].push_back(y ^ 1), adj[x ^ 1].push_back(y);
adj[y].push_back(x ^ 1), adj[y ^ 1].push_back(x);
}
bool solve() {
for (int i{0}; i < n; ++i) {
if (mark[2 * i] || mark[2 * i + 1]) continue;
while (!s.empty()) s.pop();
if (!dfs(2 * i)) {
while (!s.empty()) mark[s.top()] = false, s.pop();
if (!dfs(2 * i + 1)) return false;
}
}
return true;
}
};
```
:::
### Eulerian cycle (TODO)
```cpp!
vector<int> euler_cycle(vector<vector<pair<int, int>>>& adj, int w = 0) {
int n{adj.size()}, m{};
for (int v{0}; v < n; ++v) m += adj[v].size();
m /= 2;
vector<int> res{};
stack<int> stk{};
stk.push(w);
vector<int> nxt(n);
vector<bool> usd(m);
while (!stk.empty()) {
auto u{stk.top()};
while (nxt[u] < adj[u].size() && usd[adj[u][nxt[u]].se]) ++nxt[u];
if (nxt[u] < adj[u].size()) {
auto [v, i]{adj[u][nxt[u]]};
++nxt[u], usd[i] = true, stk.push(v);
} else res.push_back(u), stk.pop();
}
return res;
}
```
### Shortest Path
```cpp!
// Bellman-Ford algorithm
template<typename T>
optional<vector<optional<T>>> Bellman_Ford(const vector<vector<pair<int, T>>>& adj, int s) {
const auto& n{adj.size()};
vector<optional<T>> d(n, nullopt);
d[s] = 0;
queue<int> qu{}, qu2{};
vector<bool> in(n, false), in2(n, false);
qu.push(s), in[s] = true;
for (int i{0}; i < n; ++i) { // at most n-1 edges
while (!qu.empty()) {
int u{qu.front()};
qu.pop(), in[u] = false;
for (auto& [v, w] : adj[u])
if (!d[v] || d[v] > d[u].value() + w) { // relax
d[v] = d[u].value() + w;
if (!in2[v]) qu2.push(v), in2[v] = true;
}
}
qu.swap(qu2), in.swap(in2);
}
if (qu.empty()) return d;
return nullopt; // if negative cycle
}
```
```cpp!
// Dijkstra algorithm
template<typename T>
vector<optional<T>> Dijkstra(const vector<vector<pair<int, T>>>& adj, int s) {
const auto& n{adj.size()};
vector<optional<T>> d(n, nullopt);
d[s] = 0;
vector<bool> found(n, false);
using pq_t = pair<T, int>;
priority_queue<pq_t, vector<pq_t>, greater<pq_t>> pq{};
pq.emplace(0, s);
while (!pq.empty()) {
auto [_, u]{pq.top()}; pq.pop();
if (found[u]) continue;
found[u] = true;
for (auto& [v, w] : adj[u])
if (!d[v] || d[v] > d[u].value() + w) {
d[v] = d[u].value() + w;
pq.emplace(d[v].value(), v);
}
}
return d;
}
```
```cpp!
// Floyd-Warshall algorithm
template<typename T>
vector<vector<optional<T>>> Floyd_Warshall(const vector<vector<optional<T>>>& adj) {
const auto& n{adj.size()};
auto d{adj};
for (int i{0}; i < n; ++i) d[i][i] = 0;
for (int k{0}; k < n; ++k)
for (int i{0}; i < n; ++i)
for (int j{0}; j < n; ++j) {
if (!d[i][k] || !d[k][j]) continue; // no value
if (!d[i][j] || d[i][j] > d[i][k].value() + d[k][j].value())
d[i][j] = d[i][k].value() + d[k][j].value();
}
return d;
}
```
### Matching
#### maximum cardinality bipartite matching
```cpp!
class BipartiteMatching {
vector<bool> vis;
vector<int> mx{}, my{};
bool dfs(const vector<vector<int>>& adj, int x) {
for (auto& y : adj[x]) {
if (vis[y]) continue;
vis[y] = true;
if (my[y] == -1 || dfs(adj, my[y])) {
mx[x] = y, my[y] = x;
return true;
}
}
return false;
}
public:
pair<vector<int>, vector<int>> operator()(const vector<vector<int>>& adj, int ny) {
vis.resize(ny), mx.assign(adj.size(), -1), my.assign(ny, -1);
int c{0};
for (int x{0}; x < adj.size(); ++x) {
fill(vis.begin(), vis.end(), false);
if (dfs(adj, x)) ++c;
}
return {move(mx), move(my)};
}
};
```
```cpp!
// maximum cardinality bipartite matching
// Hopcroft–Karp algorithm - O(E√V)
class bipartiteGraph {
int nx, ny;
vector<int> dx, dy;
bool dfs(int y) {
for (auto& x : adjy[y]) {
if (dx[x] + 1 != dy[y]) continue;
if (mx[x] == -1 || dfs(mx[x])) {
my[y] = x, mx[x] = y;
dy[y] = -1; // augmented - remove from BFS forest
return true;
}
}
dy[y] = -1; // no augmenting path - remove from BFS forest
return false;
}
bool bfs() {
fill(dx.begin(), dx.end(), -1);
fill(dy.begin(), dy.end(), -1);
queue<int> qu{}, qu2{};
for (int x{0}; x < nx; ++x) // use all unmatched x's as roots
if (mx[x] == -1) qu.push(x), dx[x] = 0;
bool found{false};
while (!qu.empty() && !found) { // stop at the level found
while (!qu.empty()) {
auto x{qu.front()}; qu.pop();
for (auto& y : adjx[x])
if (dy[y] == -1 && my[y] != x) {
dy[y] = dx[x] + 1;
if (my[y] == -1) found = true;
else qu2.push(my[y]), dx[my[y]] = dy[y] + 1;
}
}
qu.swap(qu2);
}
return found;
}
vector<vector<int>> adjx, adjy;
public:
vector<int> mx, my;
bipartiteGraph(int _nx, int _ny) : nx{_nx}, ny{_ny}, dx(nx), dy(ny)
, adjx(nx), adjy(ny), mx(nx, -1), my(ny, -1) {}
void addEdge(int x, int y) {
adjx[x].push_back(y), adjy[y].push_back(x);
}
int Hopcroft_Karp() {
int c{0};
while (bfs()) {
for (int y{0}; y < ny; ++y)
if (my[y] == -1 && dy[y] != -1) c += dfs(y);
}
return c;
}
};
```
#### maximum weight bipartite matching
```cpp!
// Hungarian Algorithm - O(V^3)
class bipartiteGraph {
int n;
vector<long long> lx, ly;
vector<bool> vx, vy;
queue<int> qu{}; // only X's vertices
vector<int> py;
vector<long long> dy; // smallest (lx[x] + ly[y] - w[x][y])
vector<int> pdy; // & which x
void relax(int x) {
for (int y{0}; y < n; ++y)
if (lx[x] + ly[y] - w[x][y] < dy[y])
dy[y] = lx[x] + ly[y] - w[x][y], pdy[y] = x;
}
void augment(int y) {
while (y != -1) {
int x{py[y]}, yy{mx[x]};
mx[x] = y, my[y] = x;
y = yy;
}
}
bool bfs() {
while (!qu.empty()) {
int x{qu.front()}; qu.pop();
for (int y{0}; y < n; ++y) {
if (!vy[y] && lx[x] + ly[y] == w[x][y]) {
vy[y] = true, py[y] = x;
if (my[y] == -1) return augment(y), true;
int xx{my[y]};
qu.push(x), vx[xx] = true, relax(xx);
}
}
}
return false;
}
void relabel() {
long long d{numeric_limits<long long>::max()};
for (int y{0}; y < n; ++y) if (!vy[y]) d = min(d, dy[y]);
for (int x{0}; x < n; ++x) if (vx[x]) lx[x] -= d;
for (int y{0}; y < n; ++y) if (vy[y]) ly[y] += d;
for (int y{0}; y < n; ++y) if (!vy[y]) dy[y] -= d;
}
bool expand() { // expand one layer with new equality edges
for (int y{0}; y < n; ++y) {
if (!vy[y] && dy[y] == 0) {
vy[y] = true, py[y] = pdy[y];
if (my[y] == -1) return augment(y), true;
int xx{my[y]};
qu.push(xx), vx[xx] = true, relax(xx);
}
}
return false;
}
public:
vector<vector<long long>> w;
vector<int> mx, my;
bipartiteGraph(int _n) : n{_n}, lx(n), ly(n), vx(n), vy(n), py(n), dy(n), pdy(n),
w(n, vector<long long>(n, 0)), mx(n, -1), my(n, -1) {}
long long Hungarian() {
for (int i{0}; i < n; ++i) {
lx[i] = ly[i] = 0;
for (int j{0}; j < n; ++j)
lx[i] = max(lx[i], w[i][j]);
mx[i] = my[i] = -1;
}
for (int x{0}; x < n; ++x) {
fill(vx.begin(), vx.end(), false);
fill(vy.begin(), vy.end(), false);
fill(dy.begin(), dy.end(), numeric_limits<long long>::max());
qu = {}, qu.push(x), vx[x] = true, relax(x);
while (!bfs()) {
relabel();
if (expand()) break;
}
}
long long weight{0};
for (int x{0}; x < n; ++x) weight += w[x][mx[x]];
return weight;
}
};
```
#### maximum cardinality matching
```cpp!
// maximum cardinality matching
// Micali-Vazirani algorithm
// todo
```
```cpp!
// maximum matching
// Edmonds' Blossom algorithm
// O(VEα(E,V))
class graph {
int n;
vector<int> p, d, a, c1, c2;
// (alternating tree), (unvisited: -1, even: 0, odd: 1), (auxiliary for lca), (cross edge)
int label{0};
/* DSU */
vector<int> dsu_p, dsu_sz, dsu_b;
void dsu_reset() {
fill(dsu_p.begin(), dsu_p.end(), -1);
fill(dsu_sz.begin(), dsu_sz.end(), 1);
iota(dsu_b.begin(), dsu_b.end(), 0);
}
int find(int x) { // collapsing find
return dsu_p[x] == -1 ? x : dsu_p[x] = find(dsu_p[x]);
}
int base(int x) { return dsu_b[find(x)]; }
void contract(int x, int y) { // weighted union
auto rx{find(x)}, ry{find(y)};
if (rx == ry) return ;
auto b{dsu_b[rx]}; // treat x's base as new base
if (dsu_sz[rx] < dsu_sz[ry]) swap(rx, ry);
dsu_p[ry] = rx, dsu_sz[rx] += dsu_sz[ry], dsu_b[rx] = b;
}
/*******/
queue<int> qu{}; // only even vertices
void handle_one_side(int x, int y, int b) {
for (int v{base(x)}; v != b; v = base(p[v])) {
c1[v] = x, c2[v] = y, contract(b, v);
if (d[v] == 1) qu.push(v); // odd vertex -> even vertex
}
}
int lca(int u, int v) {
++label; // use next number in order not to clear a
while (true) {
if (a[u] == label) return u;
a[u] = label;
if (p[u] != -1) u = base(p[u]);
swap(u, v);
}
}
void augment(int r, int y) {
if (r == y) return ;
if (d[y] == 0) { // even vertex -> odd vertex
// check d[y] == 0 instead of d[p[y]] == 1, so (m[y], y) is in the blossom
augment(m[y], m[c1[y]]);
augment(r, m[c2[y]]);
m[c1[y]] = c2[y], m[c2[y]] = c1[y];
} else {
int x{p[y]};
augment(r, m[x]);
m[x] = y, m[y] = x;
}
}
bool bfs(int r) {
dsu_reset();
fill(d.begin(), d.end(), -1);
qu = {}, qu.push(r), d[r] = 0;
while (!qu.empty()) {
int x{qu.front()}; qu.pop();
for (auto& y : adj[x]) {
if (base(x) == base(y)) continue;
if (d[y] == -1) {
p[y] = x, d[y] = 1;
if (m[y] == -1) { // augmenting path
augment(-1, y);
return true;
} else {
p[m[y]] = y, d[m[y]] = 0;
qu.push(m[y]);
}
} else if (d[base(y)] == 0) { // blossom
int b{lca(base(x), base(y))};
handle_one_side(x, y, b);
handle_one_side(y, x, b);
}
}
}
return false;
}
public:
vector<vector<int>> adj;
vector<int> m;
explicit graph(int _n) : n{_n}, p(n, -1), d(n), a(n), c1(n), c2(n), dsu_p(n), dsu_sz(n), dsu_b(n), adj(n), m(n, -1) {}
int Edmonds() {
int c{0};
for (int v{0}; v < n; ++v)
if (m[v] == -1 && bfs(v)) ++c;
return c;
}
};
```
#### maximum weight matching (TODO)
```cpp!
// maximum weight matching
// Edmonds' Blossom algorithm
```
#### stable matching
```cpp!
// Stable Marriage Problem (stable matching)
// Gale–Shapley algorithm
class bipartiteGraph {
public:
int n;
vector<vector<int>> px, py;
vector<int> mx, my;
bipartiteGraph(int _n) : n(_n), px(n, vector<int>(n)), py(n, vector<int>(n)), mx(n, -1), my(n, -1) {}
void Gale_Shapley() {
queue<int> qu{};
vector<priority_queue<pair<int, int>>> pq(n);
for (int x{0}; x < n; ++x) {
qu.push(x);
for (int y{0}; y < n; ++y) pq[x].emplace(px[x][y], y);
}
while (!qu.empty()) {
auto x{qu.front()}; qu.pop();
int y;
do {
y = pq[x].top().second; pq[x].pop();
} while (my[y] != -1 && py[y][x] <= py[y][my[y]]);
if (my[y] != -1) mx[my[y]] = -1, qu.push(my[y]); // y prefers x to my[y]
mx[x] = y, my[y] = x;
}
}
};
```
### Flow
#### maximum flow
```cpp!
// minimum cut maximum flow
// relabel-to-front algorithm
// (an implemention of generic push-relabel algorithm)
// handle the capacity, c, of a vertex, v.
// Add a new vertex, v', and an edge v-v' with capacity c.
template<typename T> //requires signed_integral<T>
void handle_vertex_capacity(int u, T c, int& n, vector<tuple<int, int, T>>& e) {
int w{n++};
for (auto& [_u, _v, _c] : e) if (_u == u) _u = w;
e.emplace_back(u, w, c);
}
template<typename T> //requires signed_integral<T>
class flowNetwork {
int n, s, t;
vector<int> h; // height label
vector<T> e; // excess
vector<vector<int>> adj; // adjacent list
vector<vector<T>> c, r{};
vector<size_t> p; // record the position processed
list<int> lst{}; // topological sort of admissible network
void push(int u, int v) {
T f{min(e[u], r[u][v])};
r[u][v] -= f, r[v][u] += f;
e[u] -= f, e[v] += f;
}
void relabel(int u) {
int mn{numeric_limits<int>::max()};
for (auto& v : adj[u])
if (r[u][v]) mn = min(mn, h[v]);
h[u] = mn + 1;
}
void init() {
fill(h.begin(), h.end(), 0);
fill(e.begin(), e.end(), 0);
lst.clear();
for (int u{0}; u < n; ++u) if (u != s && u != t) lst.push_back(u);
fill(p.begin(), p.end(), 0);
r = c;
}
void preflow() {
h[s] = n;
for (auto& v : adj[s]) {
e[v] += r[s][v], e[s] -= r[s][v];
r[v][s] += r[s][v], r[s][v] = 0;
}
}
bool discharge(int u) {
bool flag{false};
while (e[u]) {
for (; p[u] < adj[u].size(); ++p[u]) {
int v{adj[u][p[u]]};
if (r[u][v] && h[u] == h[v] + 1) push(u, v);
if (!e[u]) break;
}
if (e[u]) relabel(u), p[u] = 0, flag = true;
}
return flag; // return if relabel or not
}
public:
flowNetwork(int _n) : n{_n}, s{0}, t{n - 1}, h(n), e(n), adj(n), c(n, vector<T>(n)), p(n) {}
void set_st(int _s, int _t) {
s = _s, t = _t;
}
void add_edge(int u, int v, T _c) {
if (!c[u][v] && !c[v][u]) adj[u].push_back(v), adj[v].push_back(u);
c[u][v] += _c;
}
T max_flow() {
assert(s != t);
init();
preflow();
auto it{lst.begin()};
while (it != lst.end()) {
if (discharge(*it)) { // relabel-to-front
lst.push_front(*it), lst.erase(it);
it = lst.begin();
}
it = next(it);
}
return e[t];
};
};
```
## String
### prefix function (LPS)
```cpp!
// longest prefix which is also suffix
template<typename T>
vector<int> prefix_func(const vector<T>& v) {
vector<int> len(v.size());
for (int i{1}; i < v.size(); ++i) {
int j{len[i - 1]};
while (j > 0 && v[i] != v[j]) j = len[j - 1];
len[i] = j + (v[i] == v[j]);
}
return len;
}
```
#### Count occurrences of prefix
```cpp!
vector<int> cnt(n + 1);
cnt[n] = 1;
for (int i{0}; i < n; ++i) ++cnt[lps[i]];
for (int i{n}; i >= 1; --i) cnt[lps[i - 1]] += cnt[i];
```
### KMP
```cpp!
template<typename T>
vector<int> KMP(const vector<T> &s, const vector<T> &p) {
vector<int> lps{prefix_func(p)}, v{};
for (int i{0}, j{0}; i < s.size(); ++i) {
while (j > 0 && s[i] != p[j]) j = lps[j - 1];
if (s[i] == p[j] && ++j == p.size()) v.push_back(i), j = lps[j - 1];
}
return v;
}
```
#### Building automaton
```cpp!
vector<array<int, 26>> fsm(p.size() + 1);
fsm[0][p[0] - 'a'] = 1;
for (int j{1}; j <= p.size(); ++j)
for (int c{0}; c < 26; ++c) {
if (j < m && p[j] == 'a' + c) fsm[j][c] = j + 1;
else fsm[j][c] = fsm[lps[j - 1]][c];
}
```
### Z function
```cpp!
template<typename T>
vector<int> z_func(const vector<T>& v) {
int n{v.size()};
vector<int> z(n);
int l{0}, r{0};
for(int i{1}; i < n; i++) {
if(i < r) z[i] = min(z[i - l], r - i);
while (i + z[i] < n && v[z[i]] == v[i + z[i]]) ++z[i];
if (r < i + z[i]) l = i, r = i + z[i];
}
return z;
}
```
#### Count occurrences of prefix
```cpp!
z[0] = n;
vector<int> cnt(n + 1);
for (int i{0}; i < n; ++i) ++cnt[z[i]];
for (int i{n - 1}; i >= 0; --i) cnt[i] += cnt[i + 1];
```
### Suffix Array
```cpp!
template<typename T>
vector<int> suffix_array(const vector<T>& _v, T s = T{}) {
int n{_v.size()};
vector<pair<T, int>> v(n);
for (int i{0}; i < n; ++i) v[i] = {_v[i], i};
v.emplace_back(s, n++);
sort(v.begin(), v.end());
vector<int> c(n), p(n);
c[v[0].se] = 0;
for (int i{1}; i < n; ++i) c[v[i].se] = c[v[i - 1].se] + (v[i - 1].fi < v[i].fi);
for (int i{0}; i < n; ++i) p[i] = v[i].se;
for (int k{0}; (1 << k) < n; ++k) {
vector<int> t(n);
for (int i{0}; i < n; ++i) t[i] = (p[i] - (1 << k) + n) % n;
vector<int> cnt(n, 0);
for (int i{0}; i < n; ++i) ++cnt[c[t[i]]];
for (int i{1}; i < n; ++i) cnt[i] += cnt[i - 1];
for (int i{n - 1}; i >= 0; --i) p[--cnt[c[t[i]]]] = t[i];
vector<int> nc(n);
nc[p[0]] = 0;
for (int i{1}; i < n; ++i) {
nc[p[i]] = nc[p[i - 1]] + (c[p[i - 1]] < c[p[i]] ||
c[(p[i - 1] + (1 << k)) % n] < c[(p[i] + (1 << k)) % n]);
}
c.swap(nc);
}
return p.erase(p.begin()), p;
}
```
### LCP array
```cpp!
// longest common prefix of consecutive suffixes
template<typename T>
vector<int> lcp_array(const vector<T>& v, const vector<int>& p) {
int n{v.size()};
vector<int> r(n);
for (int i{0}; i < n; ++i) r[p[i]] = i;
vector<int> lcp(n - 1, 0);
for (int i{0}, k{0}; i < n; ++i, k = k ? k - 1 : 0) {
if (r[i] == n - 1) continue; // k should be 0
int j{p[r[i] + 1]};
while (i + k < n && j + k < n && v[i + k] == v[j + k]) ++k;
lcp[r[i]] = k;
}
return lcp;
}
```
### Aho–Corasick algorithm
```cpp!
template<size_t k = 26>
class AC_automaton {
struct vertex_t {
int pa;
bool exit{false};
int tag;
array<int, k> next{}, go;
int slink, elink, ecnt;
explicit vertex_t(int _pa = 0) : pa{_pa} { }
};
vector<vertex_t> t;
public:
explicit AC_automaton(const vector<pair<vector<int>, int>>& dict) : t(1) {
for (auto& [s, x] : dict) {
int v{0};
for (auto& c : s) {
if (!t[v].next[c]) {
t[v].next[c] = t.size();
t.emplace_back(v);
}
v = t[v].next[c];
}
if (v) t[v].exit = true, t[v].tag = x;
}
queue<int> qu{};
qu.push(0), t[0].slink = 0;
for (size_t c{0}; c < k; ++c) t[0].go[c] = 0;
t[0].elink = 0, t[0].ecnt = 0;
while (!qu.empty()) {
auto u{qu.front()}; qu.pop();
if (u) {
int w{t[u].slink};
for (size_t c{0}; c < k; ++c) t[u].go[c] = (t[w].next[c] ? t[w].next[c] : t[w].go[c]);
t[u].elink = (t[w].exit ? w : t[w].elink), t[u].ecnt = t[w].ecnt + t[w].exit;
}
for (size_t c{0}; c < k; ++c) {
auto v{t[u].next[c]};
if (v) qu.push(v), t[v].slink = t[u].go[c];
}
}
}
size_t size() const { return t.size(); }
bool exit(int v) const { return t[v].exit; }
int tag(int v) const { return t[v].tag; }
int next(int v, int c) const { return t[v].next[c]; }
int go(int v, int c) const { return t[v].go[c]; }
int slink(int v) const { return t[v].slink; }
int elink(int v) const { return t[v].elink; }
int ecnt(int v) const { return t[v].ecnt; }
};
```
### Manacher's algorithm
```cpp!
template<typename T>
pair<vector<int>, vector<int>> Manacher(const vector<T>& _v, T s = T{}) {
int n{_v.size()};
vector<T> v(2 * n + 1, s);
for (int i{0}; i < n; ++i) v[2 * i + 1] = _v[i];
vector<int> r(2 * n + 1, 1);
int c{0}, b{1};
for (int i{1}; i < 2 * n; ++i) {
if (i < b) r[i] = min(b - i, r[c - (i - c)]);
while (i - r[i] >= 0 && i + r[i] < 2 * n + 1
&& v[i - r[i]] == v[i + r[i]]) ++r[i];
if (i + r[i] > b) c = i, b = i + r[i];
}
vector<int> ro{}, re{};
for (int i{1}; i <= 2 * n - 1; i += 2) ro.push_back(r[i] - 1);
for (int i{2}; i <= 2 * n - 2; i += 2) re.push_back(r[i] - 1);
return {ro, re};
}
```
### SAM (Suffix Automaton)
```cpp!
class SAM {
struct vertex_t {
int len, fp, link;
bool cloned;
map<int, int> next{};
explicit vertex_t(int l, int f, int lk = -1, bool c = false, const map<int, int>& n = {}) : len{l}, fp{f}, link{lk}, cloned(c), next{n} {}
};
vector<vertex_t> st{vertex_t{0, 0}};
int lst{0};
public:
explicit SAM() {}
void add_back(int ch) {
int cur{st.size()}; st.emplace_back(st[lst].len + 1, st[lst].len);
int p{lst};
while (p != -1 && !st[p].next.count(ch)) st[p].next[ch] = cur, p = st[p].link;
if (p != -1) {
int q{st[p].next[ch]};
if (st[p].len + 1 != st[q].len) {
int c{st.size()}; st.emplace_back(st[p].len + 1, st[q].fp, st[q].link, true, st[q].next);
while (p != -1 && st[p].next[ch] == q) st[p].next[ch] = c, p = st[p].link;
st[q].link = st[cur].link = c;
} else st[cur].link = q;
} else st[cur].link = 0;
lst = cur;
}
size_t size() const { return st.size(); }
int last() const { return lst; }
int len(int v) const { return st[v].len; }
int fp(int v) const { return st[v].fp; }
int link(int v) const { return st[v].link; }
int cloned(int v) const { return st[v].cloned; }
const map<int, int>& next(int v) const { return st[v].next; }
int next(int v, int ch) const {
auto it{st[v].next.find(ch)}; return it == st[v].next.end() ? -1 : it->second;
}
};
```
#### Number of different substrings
```cpp!
long long c{};
for(int i{1}; i < sam.size(); i++) c += sam.len(i) - sam.len(sam.link(i));
```
#### Number of occurrences
```cpp!
vector<pair<int, int>> tmp(sam.size());
for (int v{0}; v < sam.size(); ++v) tmp[v] = {sam.len(v), v};
sort(tmp.rbegin(), tmp.rend());
vector<long long> cnt(sam.size());
for (int v{0}; v < sam.size(); ++v) {
if (!sam.cloned(v)) cnt[v] = 1;
}
for (auto& [len, v] : tmp) {
if (v) cnt[sam.link(v)] += cnt[v];
}
```
# Math
## General
### Integer Square Root
```cpp!
unsigned long long Intsqrt(unsigned long long x) {
unsigned long long l{0ULL}, r{1ULL << 32};
while (r - l > 1) {
auto m{(l + r) / 2};
if (m * m <= x) l = m;
else r = m;
}
return l;
}
```
### Combination
```cpp!
ll Cnm(ll n, ll m) {
if (m > n / 2) m = n - m;
ll r{1};
for (ll i{1}, j{n}; i <= m; ++i, --j) r *= j, r /= i;
return r;
}
```
### Mint
```cpp!
using ll = long long;
const ll PM{1000000007};
ll MN(ll a) { return (a % PM + PM) % PM; }
ll MA(ll a, ll b) { return (a + b) % PM; }
ll MS(ll a, ll b) { return (a - b + PM) % PM; }
ll MM(ll a, ll b) { return (a * b) % PM; }
ll MP(ll a, ll b) {
ll res{1};
do {
if (b & 1) res = MM(res, a);
} while (a = MM(a, a), b >>= 1);
return res;
}
ll MI(ll a) { return MP(a, PM - 2); }
ll MD(ll a, ll b) { return MM(a, MI(b)); }
```
:::spoiler class version
```cpp!
class mint {
long long n{0};
public:
constexpr static long long mod{set value first};
mint()=default;
mint(long long _n) : n{(_n % mod + mod) % mod} {}
mint(const mint&)=default;
mint& operator=(const mint&)=default;
operator long long() { return n; }
mint& operator+=(const mint &x) {
n = (n + x.n) % mod;
return *this;
}
mint& operator*=(const mint &x) {
n = (n * x.n) % mod;
return *this;
}
mint& operator-=(const mint &x) { return (*this) += mod - x.n; }
};
template<class L, class R>
enable_if_t<(is_integral_v<L> || is_same_v<L, mint>) &&
(is_integral_v<R> || is_same_v<R, mint>), mint>
operator+(const L& lhs, const R& rhs) {
return mint{lhs} += mint{rhs};
}
template<class L, class R>
enable_if_t<(is_integral_v<L> || is_same_v<L, mint>) &&
(is_integral_v<R> || is_same_v<R, mint>), mint>
operator-(const L& lhs, const R& rhs) {
return mint{lhs} -= mint{rhs};
}
template<class L, class R>
enable_if_t<(is_integral_v<L> || is_same_v<L, mint>) &&
(is_integral_v<R> || is_same_v<R, mint>), mint>
operator*(const L& lhs, const R& rhs) {
return mint{lhs} *= mint{rhs};
}
```
:::
### Euler's totient function
```cpp!
int euler_phi(int n) {
int res{n};
for (int i{2}; i * i <= n; ++i) {
if (n % i) continue;
while (n % i == 0) n /= i;
res = res / i * (i - 1);
}
if (n > 1) res = res / n * (n - 1);
return res;
}
```
```cpp!
constexpr int maxn{100000};
vector<int> phi{[] {
vector<int> v(maxn + 1); v[1] = 1;
for (int i{2}; i <= maxn; ++i) {
if (v[i]) continue;
v[i] = i;
for (int j{i}; j <= maxn; j += i) {
if (!v[j]) v[j] = j;
v[j] = v[j] / i * (i - 1);
}
}
return v;
}()};
```
### Extended Euclidean Algorithm
:::spoiler original version
> [Explanation](https://hackmd.io/Whu1Q5dcSNuMMHXLwM5yzg)
```cpp!
// Extended Eculidean Algorithm
// solve Linear Diophantine Equation
template<typename T>
pair<T, T> operator- (const pair<T, T>& x, const pair<T, T>& y) {
return {x.first - y.first, x.second - y.second};
}
template<typename T>
pair<T, T> operator* (const T& k, const pair<T, T>& p) {
return {k * p.first, k * p.second};
}
template<typename T>
pair<T, T> extEuc(T x, T y) {
pair<T, T> p{1, 0}, q{0, 1};
if (x < y) swap(x, y), swap(p, q); // let x >= y
while (y) {
auto n{p - (x / y) * q};
auto r{x % y};
x = y, p = q;
y = r, q = n;
}
return p;
}
```
:::
```cpp!
/* solve x, y for ax + by = gcd(a, b) = g */
template<typename T>
void extEcu(T a, T b, T &g, T &x, T &y) {
if (b) extEcu(b, a % b, g, y, x), y -= x * (a / b);
else g = a, x = 1, y = 0;
}
```
### Fast Power
```cpp!
/* Iterative Function to calculate (a^b) % mod in O(log b) */
long long power_mod(long long a, long long b, long long mod) {
long long res{1};
while (b > 0) {
if (b & 1) res = res * a % mod;
b >>= 1;
a = (a * a) % mod;
}
return res;
}
```
### Modular Multiplicative Inverse
```cpp!
/* exists when a and mod are coprime */
long long MI(long long a, long long mod) {
return power_mod(a, euler_phi(mod) - 1, mod);
}
```
```cpp!
long long MI(long long a, long long mod) {
long long d, x, y;
extEcu(a, mod, d, x, y);
return d == 1 ? (x + mod) % mod : -1;
}
```
### Chinese Remainder Theorem
```cpp!
// coprime (p^k)
pair<ll, ll> CRT(const vector<pair<ll, ll>>& congruences) {
ll M{1}, sol{};
for (auto& [m, a] : congruences) M *= m;
for (auto& [m, a] : congruences) {
ll x{M / m}, y{MI(x, m)};
sol = MA(sol, MM(MM(a, x, M), y, M), M);
}
return {M, sol};
}
```
### Floyd's cycle-finding algorithm
```cpp!
template<class value_t, class func_t>
pair<size_t, size_t> FloydsTortoiseAndHare(value_t x0, func_t func) {
value_t t{x0}, h{x0};
do { t = func(t), h = func(func(h)); } while (t != h);
size_t mu{0};
t = x0;
while (t != h) t = func(t), h = func(h), ++mu;
size_t lam{0};
do { h = func(h), ++lam; } while (t != h);
return {mu, lam};
}
```
## Computational Geometry
```cpp!
struct Point {
double x, y;
Point(double _x = 0, double _y = 0) : x{_x}, y{_y} {}
};
using Vector = Point;
Vector operator+(const Vector &v1, const Vector &v2) { return {v1.x + v2.x, v1.y + v2.y}; }
Vector operator-(const Vector &v1, const Vector &v2) { return {v1.x - v2.x, v1.y - v2.y}; }
Vector operator*(const Vector &v, double d) { return {v.x * d, v.y * d}; }
Vector operator/(const Vector &v, double d) { return {v.x / d, v.y / d}; }
bool operator<(const Point &p1, const Point &p2) {
if (p1.x != p2.x) return p1.x < p2.x;
return p1.y < p2.y;
}
int dcmp(double x) {
if (fabs(x) < 1e-12) return 0;
return x > 0 ? 1 : -1;
}
bool operator==(const Point &p1, const Point &p2) {
return dcmp(p1.x - p2.x) == 0 && dcmp(p1.y - p2.y) == 0;
}
```
### Basic
```cpp!
double dot(const Vector &v1, const Vector &v2) { return v1.x * v2.x + v1.y * v2.y; }
double length(const Vector &v) { return sqrt(dot(v, v)); }
double cross(const Vector &v1, const Vector &v2) { return v1.x * v2.y - v1.y * v2.x; }
double angle(const Vector &v1, const Vector &v2) {
return acos(dot(v1, v2) / length(v1) / length(v2));
}
// 2 * area of the triangle {p1, p2, p3}
double area2(const Point &p1, const Point &p2, const Point &p3) {
return cross(p2 - p1, p3 - p1);
}
Vector rotate(const Vector &v, double rad) {
return {v.x * cos(rad) - v.y * sin(rad),
v.x * sin(rad) + v.y * cos(rad)};
}
```
### Intersection of Two Lines
```cpp!
// get the intersection of 2 lines below:
// 1. p + t1 * v
// 2. q + t2 * u
Point LineIntersection(const Point &p, const Vector &v,
const Point &q, const Vector &u) {
Vector w = q - p;
double t1 = cross(u, w) / cross(u, v);
return p + v * t1;
}
```
### Distance to Line
```cpp!
double DistanceToLine(const Point &p, const Point &a, const Point &b) {
return abs(cross(p - a, b - a) / length(b - a));
}
double DistanceToSegment(const Point &p, const Point &a, const Point &b) {
Vector v1 = p - a, v2 = b - a, v3 = p - b;
if (a == b) return length(v1);
if (dot(v1, v2) < 0) return length(v1);
if (dot(v2, v3) > 0) return length(v3);
return DistanceToLine(p, a, b);
}
```
### Projection
```cpp!
Point Projecttion(const Point &p, const Point &a, const Point &b) {
Vector v = b - a;
return a + v * (dot(v, p - a) / dot(v, v));
}
```
### Segment Intersect
```cpp!
template<typename T>
pair<T, T> operator-(const pair<T, T>& lhs, const pair<T, T>& rhs) {
auto& [lf, ls]{lhs}; auto& [rf, rs]{rhs};
return {lf - rf, ls - rs};
}
template<typename T>
int rotation(const pair<T, T>& o, const pair<T, T>& p1, const pair<T, T>& p2) {
auto v1{p1 - o}, v2{p2 - o};
// 1 for CW, -1 for CCW, and 0 for parallel
auto& [x1, y1]{v1}; auto& [x2, y2]{v2};
auto cp{x1 * y2 - x2 * y1}; // cross product
return cp ? cp / abs(cp) : 0;
}
template<typename T>
bool intersect(pair<pair<T, T>, pair<T, T>> s1, pair<pair<T, T>, pair<T, T>> s2) {
auto& [p1, p2]{s1}; auto& [p3, p4]{s2};
auto r1{rotation(p1, p2, p3)};
auto r2{rotation(p1, p2, p4)};
auto r3{rotation(p3, p4, p1)};
auto r4{rotation(p3, p4, p2)};
// assume p is on line e1-e2, check if p is on segment e1-e2
auto on_segment = [](const pair<T, T>& e1, const pair<T, T>& e2, const pair<T, T>& p) {
auto& [x1, y1]{e1}; auto& [x2, y2]{e2};
auto& [x, y]{p};
return min(x1, x2) <= x && x <= max(x1, x2) &&
min(y1, y2) <= y && y <= max(y1, y2);
};
if ((r1 && r1 == -1 * r2) && (r3 && r3 == -1 * r4)) return true;
else if (r1 == 0 && on_segment(p1, p2, p3)) return true;
else if (r2 == 0 && on_segment(p1, p2, p4)) return true;
else if (r3 == 0 && on_segment(p3, p4, p1)) return true;
else if (r4 == 0 && on_segment(p3, p4, p2)) return true;
return false;
}
```
### Polygon Area
```cpp!
double PolygonArea(const vector<Point> &p) {
int n{p.size()};
double area{0};
for (int i{2}; i < n; ++i) area += cross(p[i - 1] - p[0], p[i] - p[0]);
return fabs(area / 2);
}
```
### Circles
```cpp!
struct Circle {
Point c;
double r;
Point getPoint(double a) const {
return {c.x + cos(a) * r, c.y + sin(a) * r};
}
};
int GetLineCircleIntersection(const Point& p, const Vector& v, const Circle& c, double& t1, double& t2, vector<Point>& sol) {
double c1{v.x}, c2{p.x - c.c.x}, c3{v.y}, c4{p.y - c.c.y};
double c5{c1 * c1 + c3 * c3}, c6{2 * (c1 * c2 + c3 * c4)}, c7{c2 * c2 + c4 * c4 - c.r * c.r};
double delta{c6 * c6 - 4 * c5 * c7};
if (dcmp(delta, 0) < 0) return 0;
if (dcmp(delta, 0) == 0) {
t1 = t2 = -c6 / (2 * c5), sol.push_back(p + v * t1);
return 1;
}
t1 = (-c6 - sqrt(delta)) / (2 * c5), sol.push_back(p + v * t1);
t2 = (-c6 + sqrt(delta)) / (2 * c5), sol.push_back(p + v * t2);
return 2;
}
int GetCircleIntersection(const Circle& c1, const Circle& c2, vector<Point>& sol) {
double d{Length(c1.c - c2.c)};
if (dcmp(d, 0) == 0) {
if (dcmp(c1.r, c2.r) == 0) return -1;
return 0;
}
if (dcmp(c1.r + c2.r, d) < 0) return 0;
if (dcmp(abs(c1.r - c2.r), d) > 0) return 0;
double a{Angle(c1.c - c2.c)};
double da{acos((c1.r * c1.r + d * d - c2.r * c2.r) / (2 * c1.r * d))};
Point p1{c1.getPoint(a - da)}, p2{c1.getPoint(a + da)};
sol.push_back(p1);
if (p1 == p2) return 1;
sol.push_back(p2);
return 2;
}
int getTangents(const Point& p, const Circle& c, vector<Vector>& sol) {
Vector u{c.c - p};
double d{Length(u)};
if (dcmp(d, c.r) > 0) return 0;
else if (dcmp(d, c.r) == 0) {
sol.push_back(Rotate(u, PI / 2));
return 1;
}
}
```
## Game Theory
:::info
sg value: **mex** of successive states
nim value: **xor-sum** of all independent games' sg values
:::
## Formula
### Power Summation
```cpp!
// 1 + 4 + ... + n^2
n * (n + 1) * (2 * n + 1) / 6
```
```cpp!
// 4 + 16 + ... + (2n)^2
(2 * n) * (n + 1) * (2 * n + 1) / 3
```
```cpp!
// 1 + 9 + ... + (2n+1)^2
n * (2 * n - 1) * (2 * n + 1) / 3
```
```cpp!
// 1 + 8 + ... + n^3
(n^4 + 2 * n^3 + n^2) / 4
```
```cpp!
// 1 + 16 + ... + n^4
(6 * n^5 + 15 * n^4 + 10 * n^3 - n) / 30
```
```cpp!
// 1 + 32 + ... + n^5
(2 * n^6 + 6 * n^5 + 5 * n^4 - n^2) / 12
```
```cpp!
// 1 + 64 + ... + n^6
(6 * n^7 + 21 * n^6 + 21 * n^5 - 7 * n^3 + n) / 42
```
```cpp!
// 1 + 128 + ... + n^7
(3 * n^8 + 12 * n^7 + 14 * n^6 - 7 * n^4 + 2 * n^2) / 24
```
```cpp!
// 1 + + ... + n^8
(10 * n^9 + 45 * n^8 + 60 * n^7 - 42 * n^5 + 20 * n^3 - 3 * n) / 90
```
```cpp!
// 1 + + ... + n^9
(2 * n^10 + 10 * n^9 + 15 * n^8 - 14 * n^6 + 10 * n^4 - 3 * n^2) / 20
```
```cpp!
// 1 + + ... + n^10
(6 * n^11 + 33 * n^10 + 55 * n^9 - 66 * n^7 + 66 * n^5 - 33 * n^3 + 5 * n) / 66
```
### Summation of Combination
```cpp
// Odd term
// k = 1 - n % 2
// C(1, n) + C(3, n) + ... + C(n - k, n) = 2^(n - 1)
```
# Resources
## Online Judge
* [CSES](https://cses.fi/problemset/list/)
* [UVa](https://onlinejudge.org/index.php)
* [程式自學平台](https://e-tutor.itsa.org.tw/e-Tutor/Question_bank.php?id=30)
## Reference
* [Competitive Programmer’s Handbook](https://cses.fi/book/book.pdf?fbclid=IwAR0Gru8WZO9AnOIAiPXuL9Kdk-s__QuYvSfsfgPMfjtkFLznItnyIe91wuE)
* [演算法筆記](http://web.ntnu.edu.tw/~algo/)
* [CP-Algorithms](https://cp-algorithms.com/)
* [koosaga's github](https://github.com/koosaga/olympiad/tree/master/Library)
### C++
* [cppreference](https://en.cppreference.com/w/)
* [cplusplus](http://www.cplusplus.com/)
###### tags: `cp`
Sign in with Wallet
Connect another wallet