# 競程魔法書
## C++模板
**vimrc**
``` vim=
set encoding=utf-8
set nu
syntax on
set ai
set shiftwidth=4
set tabstop=4
set expandtab
set smartindent
set backspace=2
set ruler
set ic
set number
set mouse=a
set hlsearch
inoremap {<CR> {<CR>}<Esc>ko
inoremap " ""<Esc>i
inoremap ' ''<Esc>i
```
**模板**
```cpp=
#include <bits/stdc++.h>
using namespace std;
#define INITIO() ios_base::sync_with_stdio(false);cin.tie(NULL);cout<<setprecision(numeric_limits<double>::max_digits10+1);
#define MEM(i,j) memset(i,j,sizeof i);
#define X first
#define Y second
#define pb push_back
#define mp make_pair
#define ALL(a) a.begin(),a.end()
typedef long long ll;
typedef pair<int ,int> pii;
const double Eps = 1e-8;
const ll Inf = 1e18+7;
const int Mod = 1e9+7;
int main (){
INITIO();
int tt;
cin >> tt;
while (tt--) {
}
return 0;
}
```
## STL
### Map STL
```cpp=
map<T, T> m;
unordered_map<T, T> umap; //c++11
//放東西
m.insert();
m[] = ;
//找東西
auto iter = m.find();
if(iter != m.end())
cout<<iter->second<<endl;
```
### Queue STL
```cpp=
queue<T> q;
priority_queue<T> pq;
q.push();
q.pop();
q.front();
q.back();
```
### 排列組合 C(n,r)
>stackoverflow 來的，所以好像有overflow可能
```cpp=
long long C(int n, int r) {
if(r > n - r) r = n - r;
long long ans = 1;
int i;
for(i = 1; i <= r; i++) {
ans *= n - r + i;
ans /= i;
}
return ans;
}
```
### Sort & Cmp
>**sort** 預設升冪 , 改降冪使用 **greater<T>()**
```cpp=
bool myfunction (int a,int b) {
return (a<b);
}
sort(All(s), myfunction);
//用lambda，A是某種struct
sort( values.begin( ), values.end( ), [ ]( const A& l, const A& r )
{
return l.key < r.key;
});
```
## 快速冪
>一般的
```cpp=
ll Power(ll n,ll time){
ll a = 1;
while(time){
if(time & 1) a = a*n % Mod;
n = n*n % Mod;
time>>=1;
}
return a;
}
```
### 矩陣快速冪
>費式數列 , 比**dp**解更快
```cpp=
ll Pow(ll t){
t--;
vector<ll> a = {1,1,1,0},b = {1,1,1,0};
while(t){
vector<ll> c = a,d = b;
if(t&(ll)1){
for(int i=0;i<2;i++)
for(int j=0;j<2;j++)
a[2*i+j] = c[2*i]*b[j]+c[2*i+1]*b[2+j];
}
for(int i=0;i<2;i++)
for(int j=0;j<2;j++)
b[2*i+j] = d[2*i]*d[j]+d[2*i+1]*d[2+j];
t>>=(ll)1;
}
return a[0];
}
```
## 質數
>線性篩法
```cpp=
ll n;
cin >> n;
bool arr[n+1];
MEM(arr,true);
for(int i=2;i<=sqrt(n);++i)
if(arr[i])
for(int j = i*i ; j<=n ;j+=i)
arr[j] = false;
for(int k=1;k<=n;k++)
if(arr[k]) cout << k <<endl;
```
## 最小生成樹 MST
>求長度應用例
>Kruskal演算法 , 路徑壓縮沒啟發式合併
```cpp=
const int N = 5e3+5,M = 2e5+5;
int boss[N]={0};
struct Edge{
int a;
int b;
int val;
} e[M];
bool cmp(Edge a,Edge b){
return a.val<b.val;
}
int find(int idx){
if(boss[idx]==idx) return idx;
return boss[idx] = find(boss[idx]);
}
```
``` cpp=
int n,m,length=0,en=0;
cin>>n>>m;
for(int i=1;i<n+1;i++){
boss[i]=i;
}
for(auto &it : e){
cin>>it.a>>it.b>>it.val;
}
sort(e,e+m,cmp);
for(int j=0;j<m;j++){
if(find(e[j].a) == find(e[j].b)) continue;
boss[find(e[j].a)] = find(e[j].b);
length+=e[j].val;
en++;
}
if(en<n-1) cout << "orz"<<endl;
else cout <<length<<endl;
```
### Prim’s algorithm
```cpp=
#define MAX 10000
struct node {
int id,cost;
};
vector<node> E[MAX];
struct edge {
int u,v,cost;
bool operator<(const edge& y) const { return cost > y.cost; }
};
priority_queue<edge> Q;
bool vis[MAX];
//use
{
// init
for(int i = 0; i < n; i++)
vis[i] = false;
while(m--) {
cin >> u >> v >> cost;
E[u].push_back({v,cost});
E[v].push_back({u,cost});
}
// Prim's Algorithm
Q.push({0,0,0});
while(!Q.empty()) {
edge e = Q.top();
Q.pop();
u = e.v;
if(!vis[u]) {
ans += e.cost;
for(auto v : E[u]) {
if(!vis[v.id])
Q.push({u,v.id,v.cost});
}
}
vis[u] = true;
}
cout << ans << '\n';
}
```
## 最短路
**Floyd**
```cpp=
ll dp[Maxn][Maxn];
ll n,k,q;
cin >> n >> k >> q;
for(int i = 0;i<n ;i++){
for(int j =0;j<n ;j++){
if(i == j) dp[i][j] = 0;
else dp[i][j] = Inf;
}
}
while(k--){
ll a,b,w;
cin >> a >> b >> w;
dp[a][b] = w;
}
for(int i =0;i<n;i++)
for(int k=0;k<n;k++)
for(int j = k;j <n;j++)
dp[i][j] = min(dp[i][j],dp[i][k]+dp[k][j]);
```
## Dijkstra
單點搜尋
$O((E+N)log(N))$
```
struct Node{
ll n,w;
bool operator<(const node &lhs) const {
return lhs.w < w; //min heap
}
};
```
初始:
```
MEM(dist,Inf); //不連通設Inf
MEM(vis,false); //沒有拜訪
priority_queue<Node> q;
q.push((Node){1,dist[1]=0}); //看起點更改
```
循環:
```
while(q.size()){
node now = q.top();q.pop();
vis[now.n] = true;
for(auto nxt : e[now.n]){
if(vis[nxt.n]) continue;
const int &updater = dist[now.n] + nxt.w;
if(update < dist[nxt.n]) {
dist[nxt.n] = update; //鬆弛
q.push((Node){nxt.n,dist[nxt.n]});
}
}
}
```
## Bellman-Ford+檢查負環
單點搜尋
$O(NE)$
```
int w[N][N];
int d[N];
int parent[N];
void bellman_ford(int source){
for(int i=0;i<n;i++) d[i]=Inf;
d[source] = 0;
parent[source] = source;
for(int i=0;i<n-1;i++)
for(int a=0;a<n;++a)
for(int b=0;b<n;b++)
if(d[a] != Inf && w[a][b] != Inf)
if(d[a]+w[a][b] < d[b]){
d[b] = d[a] + w[a][b];
parent[b] = a;
}
}
```
負環檢查
* 如果還能鬆弛代表有負環
```
bool is_negative(){
for(int a=0;a<n;++a)
for(int b=0;b<n;++b)
if(d[a] != Inf && w[a][b] != Inf)
if(d[a]+w[a][b] < d[b])
return true;
return false;
}
```
遞迴版本 視情況修改 (一開始刻我以為是 dijkstra,但能處理處理負環)
```
ll bell(int start,int end,int val){
ll ans = Inf;
if(dis[start]<=val) return Inf;
else dis[start] = val;
if(start == end) return val;
for(int i=0;i<E[start].size();i++)
ans = min(ans,bell(E[start][i].X,end,E[start][i].Y+val));
return ans;
}
```
## Floyd-Warshell + 檢查負環
$O(N^3)$
多源搜尋
```
for(int i = 0;i<n ;i++){
for(int j =0;j<n ;j++){
if(i == j) dp[i][j] = 0;
else dp[i][j] = Inf;
}
}
...
for(int k =0;k<n;k++)
for(int i=0;i<n;i++)
for(int j = 0;j <n;j++)
dp[i][j] = min(dp[i][j],dp[i][k]+dp[k][j]);
```
負環
```
bool is_negative(){
for (int i=0; i<n; i++)
if (d[i][i] < 0) //自己到自己必須是0否則有負環
return true;
return false;
}
```
## 關節點
```cpp=
#define MAX 105
int low[MAX],dep[MAX];
bool ap[MAX];
int d = 0;
vector<int> adj[MAX];
void dfs(int u, int w == -1) {
int child{0};
low[u] = dep[u] = ++d;
for (auto& v : adj[u]) {
if (v == w) continue; // 剛剛過來的點
if (!dep[v]) { // 尚未拜訪的點
++child;
dfs(v, u);
low[u] = min(low[u], low[v]); // 更新 back edge (u 先往 child 走可能會有更好的 back edge)
if (low[v] >= dep[u]) // v 走 back edge 無法到 u 的祖先
ap[u] = true;
}
else // 已經拜訪過得的點
low[u] = min(low[u], dep[v]); // 更新 back edge (u 可能可以直接走到祖先)
}
if (w == -1)
ap[u] = (child > 1); // 檢查 root 是否為關節點 (子樹>1)
--d;
}
```
## 最低共同祖先
```cpp=
// pointer jumping technique (binary lifting)
#define MAX 30001
vector<int> node[MAX];
int DPWay[MAX][15]; /*點往上走 2^i level 所到的點，2^15=32768 > 30001*/
int treeLevel[MAX]; /*點高度*/
int treePath[MAX]; /*紀錄當前各level的點*/
void dfs(int now, int pre, int level) {
treeLevel[now] = level;
treePath[level] = now;
for(int i = 1; (1 << i) < level; i++)
DPWay[now][i] = treePath[level - (1 << i)];
for(int i = 0; i < node[now].size(); i++) {
if(node[now][i].first == pre)
continue;
DPWay[node[now][i].first][0] = now;
dfs(node[now][i].first, now, level + 1);
}
}
int LCA(int A,int B) {
int levelA = treeLevel[A], levelB = treeLevel[B];
for(int i = (int)log2(abs(levelA - levelB)); i >= 0; --i) {
if(abs(levelA - levelB) >= (1 << i)) {
if(levelA > levelB) {
A = DPWay[A][i];
levelA -= (1 << i);
}
else {
B = DPWay[B][i];
levelB -= (1 << i);
}
}
}
// A B 走到同一個高度 (假設樹高小於32，我們先看距離有沒有超過16，有的話就往上走16步；剩下的距離一定小於16，我們再看它有沒有超過8；以此類推，4,2,1，就到終點了。)
if (A == B)
return A;
for (int i = (int)log2(levelA); i >= 0; --i)
if (DPWay[A][i] != DPWay[B][i]) {
A = DPWay[A][i];
B = DPWay[B][i];
}
// u, v 一起走到 lca(u, v) 的下方
return DPWay[A][0];
// 回傳 lca(u, v)
}
```
## 線段樹
```cpp=
template <typename value_t>
class SGT {
public:
int n;
vector<value_t> t;
vector<optional<value_t>> lz;
int left(int tv) { return tv + 1; }
int right(int tv,int tl, int tm) { return tv + (tm-tl) * 2 ; }
void build(const vector<value_t>& v, int tv, int tl, int tr) {
if(tr - tl == 1) {
t[tv] = v[tl];
}
else {
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)]);
}
}
value_t merge(const value_t& x, const value_t& y) {
return 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] += len * x;
}
void push(int tv, int tl, int tr) {
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 modify(int l, int r, const value_t& x, int tv, int tl, int tr) {
if(l == tl && r == tr)
update(tv,tr-tl,x);
else {
push(tv,tl,tr);
int tm = (tl+tr) / 2;
if(r <= tm)
modify(l,r,x,left(tv),tl,tm);
else if(l >= tm)
modify(l,r,x,right(tv,tl,tm),tm,tr);
else {
modify(l,tm,x,left(tv),tl,tm);
modify(tm,r,x,right(tv,tl,tm),tm,tr);
}
t[tv] = merge(t[left(tv)],t[right(tv,tl,tm)]);
}
}
value_t query(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 query(l,r,left(tv),tl,tm);
else if(l >= tm)
return query(l,r,right(tv,tl,tm),tm,tr);
else
return merge(query(l,tm,left(tv),tl,tm),query(tm,r,right(tv,tl,tm),tm,tr));
}
public:
SGT(const vector<value_t>& v) : n{v.size()}, t(2*n), lz(2*n) {
build(v,1,0,n);
}
void modify(int l, int r, const value_t& x) { modify(l,r,x,1,0,n); }
value_t query(int l,int r) { return query(l,r,1,0,n); }
};
// a[0] a[1] a[2]
// 1 2 3
// l = 1 r = 3 k = 5
// query(l-1,r) modify(l-1,r,k)
```
### 懶人標記
```cpp=
//新值是value，修改區間是[l,r]，節點v的範圍是[L,R]
void modify(type value, int l, int r, int L, int R, vertex v){
if(l == L && r == R){
打懶標在v上;
return;
}
int M = (L + R) / 2;
if(r <= M) modify(value, l, r, L, M, v的左子節點);
else if(l > M) modify(value, l, r, M + 1, R, v的右子節點);
else{
modify(value, l, M, L, M, v的左子節點);
modify(value, M + 1, r, M + 1, R, v的右子節點);
}
用兩個子節點的答案更新v的答案;
}
```
修改
```cpp=
//查詢範圍是[l,r]，節點v的區間範圍是[L,R]
type query(int l, int r, int L, int R, vertex v){
if(l == L && r == R)
return ans of v;
把節點v的懶標下推; //<----
int M = (L + R) / 2;
if(r <= M)
return query(l, r, L, M, v的左子節點);
else if(l > M)
return query(l, r, M + 1, R, v的右子節點);
else{
return 用query(l, M, L, M, v的左子節點)
和query(M + 1, r, M + 1, R, v的右子節點)
得出區間[l,r]的答案;
}
}
```
## Computation Geometry
```cpp=
struct Point{
double x, y;
Point(double x = 0, double y = 0) : x(x), y(y){}
Point operator+(const Point &b)const{ // vector addition
return Point(x + b.x, y + b.y);
}
Point operator-(const Point &b)const{ // vector subtraction
return Point(x - b.x, y - b.y);
}
Point operator*(double d)const{ // scale (multiplication)
return Point(x * d, y * d);
}
Point operator/(double d)const{ // scale (division)
return Point(x / d, y / d);
}
double dot(const Point &b)const{ // vector dot
return x * b.x + y * b.y;
}
double cross(const Point &b)const{ // vector cross
return x * b.y - y * b.x;
}
};
// premise: three points are collinear
bool btw(const Point &p1, const Point &p2, const Point &p3){
return (p1 - p3).dot(p2 - p3) <= 0;
}
bool collinearity(const Point &p1, const Point &p2, const Point &p3){
return (p1 - p3).cross(p2 - p3) == 0;
}
bool pointOnSegment(const Point &p1, const Point &p2, const Point &p3){
return collinearity(p1, p2, p3) && btw(p1, p2, p3);
}
int ori(const Point &p1, const Point &p2, const Point &p3){
double d = (p2 - p1).cross(p3 - p1);
if(d == 0)
return 0;
return d > 0 ? 1 : -1;
}
bool seg_intersect(const Point &p1, const Point &p2, const Point &p3, const Point &p4){
int a123 = ori(p1, p2, p3);
int a124 = ori(p1, p2, p4);
int a341 = ori(p3, p4, p1);
int a342 = ori(p3, p4, p2);
if(a123 == 0 && a124 == 0)
return btw(p1, p2, p3) || btw(p1, p2, p4) || btw(p3, p4, p1) || btw(p3, p4, p2);
return a123 * a124 <= 0 && a341 * a342 <= 0;
}
Point intersect(const Point &p1, const Point &p2, const Point &p3, const Point &p4){
int a123 = (p2 - p1).cross(p3 - p1);
int a124 = (p2 - p1).cross(p4 - p1);
return (p4 * a123 - p3 * a124) / (a123 - a124);
}
```
## 逆序數對
```cpp=
int merge_sort(vector<int> &a,int l,int r){
if(l==r) return 0;
int mid = (l+r)/2;
int lans = merge_sort(a,l,mid);
int rans = merge_sort(a,mid+1,r);
vector<int> b,c;
for(int i = l;i<=mid;i++) b.pb(a[i]);
for(int i = mid+1;i<=r;i++) c.pb(a[i]);
int i=0,j=0,ans=0;
for(int k=l;k<=r;k++){
int Min = 0;
if(i == b.size()) Min = c[j++];
else if(j == c.size()) Min = b[i++];
else if(b[i] <= c[j]) Min = b[i++];
else {
Min = c[j++];
ans += mid - i - l + 1;
}
a[k] = Min;
}
return lans+rans+ans;
}
```
## DP
### 01背包問題
```cpp=
int v[n];
int c[n];
int dp[101]={0};
int ans=0;
for(int i=0;i<n;i++){
cin>> v[i] >> c[i];
for(int j=100;j-v[i]>=0;j--){
dp[j] = max(dp[j],dp[j-v[i]]+c[i]);
ans = max(ans,dp[j]);
}
}
cout << ans <<endl;
```
## GCD
>GCD求解
```cpp=
int gcd(int a, int b)
{ return a? gcd(b%a, a) : b;}
```
### 貝祖等式 (Bezout’s Thm)
>對整數 **a,b** , 存在整數 **x,y** 使得 **ax+by = gcd(a,b)**
### 擴展歐基里德
令 **g = gcd(0,g) = gcd(a,b)**
```cpp=
int gcd(int a, int b, int &x, int &y) {
if(!a) {
x = 0; // x 為任意整數
y = 1;
return b;
}
int g = gcd(b%a, a, x, y);
int xp = y - b/a * x, yp = x; // p := previous
x = xp, y = yp;
return g;
}
```
**exgcd** 可求模反元素 , 但僅於 b 跟 m 互質
## Karatsuba algorithm
```cpp=
// 不要爆掉，感覺都是對的
int getLength(long long value) {
int counter = 0;
while (value != 0) {
counter++;
value /= 10;
}
return counter;
}
long long karatsuba(long long x, long long y) {
int xLength = getLength(x);
int yLength = getLength(y);
// the bigger of the two lengths
int N = (int)(fmax(xLength, yLength));
// if the max length is small it's faster to just flat out multiply the two nums
if (N < 10)
return x * y;
//max length divided and rounded up
N = (N/2) + (N%2);
long long multiplier = pow(10, N);
long long b = x/multiplier;
long long a = x - (b * multiplier);
long long d = y / multiplier;
long long c = y - (d * N);
long long z0 = karatsuba(a,c);
long long z1 = karatsuba(a + b, c + d);
long long z2 = karatsuba(b, d);
return z0 + ((z1 - z0 - z2) * multiplier) + (z2 * (long long)(pow(10, 2 * N)));
}
//>
```
## 流
```cpp=
struct MaxFlow{
struct edge{
int to, cap, flow, rev;
};
vector<edge> v[MAXN];
int s, t, dis[MAXN], cur[MAXN];
int dfs(int u, int cap){
if(u == t || !cap)
return cap;
for(int &i = cur[u]; i < v[u].size(); i++){
edge &e = v[u][i];
if(dis[e.to] == dis[u] + 1 && e.flow != e.cap){
int df = dfs(e.to, min(e.cap - e.flow, cap));
if(df){
e.flow += df;
v[e.to][e.rev].flow -= df;
return df;
}
}
}
dis[u] = -1;
return 0;
}
bool bfs(){
memset(dis, -1, sizeof(dis));
queue<int> q;
q.push(s), dis[s] = 0;
while(!q.empty()){
int tmp = q.front();
q.pop();
for(auto &u : v[tmp]){
if(!~dis[u.to] && u.flow != u.cap){
q.push(u.to);
dis[u.to] = dis[tmp] + 1;
}
}
}
return dis[t] != -1;
}
int maxflow(int _s, int _t){
s = _s, t = _t;
int flow = 0, df;
while(bfs()){
memset(cur, 0, sizeof(cur));
while(df = dfs(s, INF)){
flow += df;
}
}
return flow;
}
void addedge(int st, int ed, int cap){
v[st].push_back(edge{ed, cap, 0, v[ed].size()});
v[ed].push_back(edge{st, 0, 0, v[st].size() - 1});
}
};
```
## 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) {
if (mx[x] != -1) continue;
fill(vis.begin(), vis.end(), false);
if (dfs(adj, x)) ++c;
}
return {move(mx), move(my)};
}
};
```
```cpp=
// Hopcroft–Karp algorithm - O(E√V)
class BipartiteMatching {
vector<vector<int>> adjy;
vector<int> dx, dy, mx, my;
bool dfs(const vector<vector<int>>& adjx, int y) {
for (auto& x : adjy[y]) {
if (dx[x] + 1 != dy[y]) continue;
if (mx[x] == -1 || dfs(adjx, mx[x])) {
my[y] = x, mx[x] = y;
dy[y] = -1; // remove from bfs forest
return true;
}
}
dy[y] = -1; // remove from bfs forest
return false;
}
bool bfs(const vector<vector<int>>& adjx) {
fill(dx.begin(), dx.end(), -1);
fill(dy.begin(), dy.end(), -1);
queue<int> qu{}, qu2{};
for (int x{0}; x < adjx.size(); ++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;
}
public:
pair<vector<int>, vector<int>> operator()(const vector<vector<int>>& adjx, int ny) {
adjy.assign(ny, {});
for (int x{0}; x < adjx.size(); ++x)
for (auto& y : adjx[x]) adjy[y].push_back(x);
dx.resize(adjx.size()), dy.resize(ny), mx.assign(adjx.size(), -1), my.assign(ny, -1);
int c{0};
while (bfs(adjx)) {
for (int y{0}; y < ny; ++y)
if (my[y] == -1 && dy[y] != -1) c += dfs(adjx, y);
}
return {move(mx), move(my)};
}
};
```
### 一般圖最大匹配
```cpp=
// Edmonds' Blossom Algorithm - O(VEa(V))
class Matching {
int label{0};
vector<int> p{}, d{}, a{}, c1{}, c2{}, m{};
// (alternating tree), (unvisited: -1, even: 0, odd: 1), (auxiliary for lca), (cross edge)
/* 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(const vector<vector<int>>& adj, 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<int> operator()(const vector<vector<int>>& adj) {
label = 0, p.assign(adj.size(), -1), d.resize(adj.size()), a.resize(adj.size()), c1.resize(adj.size()), c2.resize(adj.size()), m.assign(adj.size(), -1);
dsu_p.resize(adj.size()), dsu_sz.resize(adj.size()), dsu_b.resize(adj.size());
int c{0};
/*
// find random greey matching first
vector<int> X(n); iota(X.begin(), X.end(), 0);
random_shuffle(X.begin(), X.end());
for (auto& x : X) {
if (m[x] != -1) continue;
for (auto& y : adj[x])
if (m[y] == -1) {
m[x] = y, m[y] = x, ++c;
break;
}
}
*/
for (int v{0}; v < adj.size(); ++v)
if (m[v] == -1 && bfs(adj, v)) ++c;
return move(m);
}
};
```
### 二分圖最大權匹配
```cpp=
// Hungarian Algorithm - O(V^3)
class WeightedBipartiteMatching {
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
vector<int> mx{}, my{};
void relax(const vector<vector<long long>>& w, 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(const vector<vector<long long>>& w) {
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(w, 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(const vector<vector<long long>>& w) { // 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(w, xx);
}
}
return false;
}
public:
pair<vector<int>, vector<int>> operator()(const vector<vector<long long>>& w) {
n = w.size(), lx.assign(n, 0), ly.assign(n, 0), vx.resize(n), vy.resize(n), py.resize(n), dy.resize(n), pdy.resize(n), mx.assign(n, -1), my.assign(n, -1);
for (int i{0}; i < n; ++i)
for (int j{0}; j < n; ++j)
lx[i] = max(lx[i], w[i][j]);
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(w, x);
while (!bfs(w)) {
relabel();
if (expand(w)) break;
}
}
return {move(mx), move(my)};
}
};
```
## 其他
### 一些 g++ 內建的 bit function
```cpp=
__builtin_clz() //count leading 0s
__builtin_ctz() //count trailing 0s
__builtin_popcount() //count bit 1
```
### 凸包、旋轉尺
```cpp=
vector<Point> convexHull(vector<Point> v){
int sz = 0;
vector<Point> ans;
sort(v.begin(), v.end(), cmp);
for(size_t i = 0; i < v.size(); i++){
while(sz >= 2 && (ans[sz - 1] - ans[sz - 2]).cross(v[i] - ans[sz - 2]) <= 0)
ans.pop_back(), sz--;
ans.push_back(v[i]), sz++;
}
for(int i = int(v.size()) - 2, t = sz + 1; i >= 0; i--){
while(sz >= t && (ans[sz - 1] - ans[sz - 2]).cross(v[i] - ans[sz - 2]) <= 0)
ans.pop_back(), sz--;
ans.push_back(v[i]), sz++;
}
if(v.size() > 1)
ans.pop_back();
return ans;
}
double rotatingCaliper(const vector<Point> &v){
double ans = 0;
for(int i = 0; j = 0; i < v.size(); i++){
while(abs2(v[i],v[j]) <= abs2(v[i], v[(j+1)%v.size()]))
j = (j+1) % v.size();
ans = max(ans, abs2(v[i], v[j]));
}
}
```
Sign in with Wallet
Connect another wallet