Soongsil University – PS akgwi Page 1 of 25
3.23 O((V + E) log V ) Dominator Tree . . . . . . . . . . . . . . 12 7 Notes 24
3.24 O(V E) Vizing Theorem . . . . . . . . . . . . . . . . . . . . 12 7.1 Triangles/Trigonometry . . . . . . . . . . . . . . . . . . . . 24
3.25 O(E + V 3 + V 3T + V 2 2T ) Minimum Steiner Tree . . . . . 13 7.2 Calculus, Newton’s Method . . . . . . . . . . . . . . . . . . 24
Team Note of PS akgwi 3.26 O(E log V ) Directed MST . . . . . . . . . . . . . . . . . . . 13 7.3 Zeta/Mobius Transform . . . . . . . . . . . . . . . . . . . . 24
3.27 O(E log V + K log K) K Shortest Walk . . . . . . . . . . . . 13 7.4 Generating Function . . . . . . . . . . . . . . . . . . . . . . 24
3.28 O(V + E) Chordal Graph, Tree Decomposition . . . . . . . 13 7.5 Counting . . . . . . .P. . . . . . . . . . . . . . . . . . . . . 24
Jeounghui Nah, Joowon Oh, Seongseo Lee 3.29 O(V 3 ) General Matching . . . . . . . . . . . . . . . . . . . 14 7.6 Faulhaber’s Formula ( n c
k=1 k ) . . . . . . . . . . . . . . . . 25
3.30 O(V 3 ) Weighted General Matching . . . . . . . . . . . . . . 14 7.7 About Graph Degree Sequence . . . . . . . . . . . . . . . . 25
48th ICPC World Finals (Astana, Kazakhstan) 7.8 Burnside, Grundy, Pick, Hall, Simpson, Area of Quadrangle,
4 Math 15 Fermat Point, Euler, Pythagorean . . . . . . . . . . . . . . 25
4.1 Binary GCD, Extend GCD, CRT, Combination . . . . . . . 15 7.9 About Graph Minimum Cut . . . . . . . . . . . . . . . . . . 25
4.2 Partition Number . . . . . . . . . . . . . . . . . . . . . . . . 15 7.10 Matrix with Graph(Kirchhoff, Tutte, LGV) . . . . . . . . . 25
Contents 4.3 Diophantine . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 7.11 About Graph Matching(Graph with |V | ≤ 500) . . . . . . . 25
4.4 FloorSum . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 7.12 Checklist . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
1 DataStructure 1
1.1 Erasable Priority Queue . . . . . . . . . . . . . . . . . . . . 1
4.5 XOR Basis (XOR Maximization) . . . . . . . . . . . . . . . 16 1 DataStructure
4.6 Stern Brocot Tree . . . . . . . . . . . . . . . . . . . . . . . 16 1.1 Erasable Priority Queue
1.2 Convex Hull Trick (Stack, LineContainer) . . . . . . . . . . 1
1.3 Color Processor . . . . . . . . . . . . . . . . . . . . . . . . . 2 4.7 O(N 3 log 1/) Polynomial Equation . . . . . . . . . . . . . . 16
template<class T=int, class O=less<T>>
1.4 Kinetic Segment Tree . . . . . . . . . . . . . . . . . . . . . 2 4.8 Gauss Jordan Elimination . . . . . . . . . . . . . . . . . . . 16
struct pq_set {
1.5 Lazy LiChao Tree . . . . . . . . . . . . . . . . . . . . . . . . 2 4.9 Berlekamp + Kitamasa . . . . . . . . . . . . . . . . . . . . 16
priority_queue<T, vector<T>, O> q, del;
1.6 Splay Tree, Link-Cut Tree . . . . . . . . . . . . . . . . . . . 3 4.10 Linear Sieve . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
const T& top() const { return q.top(); }
4.11 Xudyh Sieve . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
int size() const { return int(q.size()-del.size()); }
2 Geometry 3 4.12 Miller Rabin + Pollard Rho . . . . . . . . . . . . . . . . . . 17
bool empty() const { return !size(); }
2.1 O(log N ) Point in Convex Polygon . . . . . . . . . . . . . . 3 4.13 Primitive Root, Discrete Log/Sqrt . . . . . . . . . . . . . . 17 void insert(const T x) { q.push(x); flush(); }
2.2 Segment Distance, Segment Reflect . . . . . . . . . . . . . . 4 4.14 Power Tower . . . . . . . . . . . . . . . . . . . . . . . . . . 17 void pop() { q.pop(); flush(); }
2.3 Tangent Series . . . . . . . . . . . . . . . . . . . . . . . . . 4 4.15 De Bruijn Sequence . . . . . . . . . . . . . . . . . . . . . . . 18 void erase(const T x) { del.push(x); flush(); }
2.4 Intersect Series . . . . . . . . . . . . . . . . . . . . . . . . . 4 4.16 Simplex / LP Duality . . . . . . . . . . . . . . . . . . . . . 18 void flush() { while(del.size() && q.top()==del.top())
2.5 O(N 2 log N ) Circles Area . . . . . . . . . . . . . . . . . . . 5 4.17 Polynomial & Convolution . . . . . . . . . . . . . . . . . . . 18 q.pop(), del.pop(); }
2.6 Segment In Polygon . . . . . . . . . . . . . . . . . . . . . . 5 4.18 Matroid Intersection . . . . . . . . . . . . . . . . . . . . . . 19 };
2.7 Polygon Cut, Center, Union . . . . . . . . . . . . . . . . . . 5
2.8 Polycon Raycast . . . . . . . . . . . . . . . . . . . . . . . . 5 5 String 20 1.2 Convex Hull Trick (Stack, LineContainer)
2.9 O(N log N ) Shamos-Hoey . . . . . . . . . . . . . . . . . . . 6 5.1 KMP, Hash, Manacher, Z . . . . . . . . . . . . . . . . . . . 20
struct Line{ // call init() before use
2.10 O(N log N ) Half Plane Intersection . . . . . . . . . . . . . . 6 5.2 Aho-Corasick . . . . . . . . . . . . . . . . . . . . . . . . . . 20 ll a, b, c; // y = ax + b, c = line index
2.11 O(M log M ) Dual Graph . . . . . . . . . . . . . . . . . . . . 6 5.3 O(N log N ) SA + LCP . . . . . . . . . . . . . . . . . . . . . 20 Line(ll a, ll b, ll c) : a(a), b(b), c(c) {}
2.12 O(N 2 log N ) Bulldozer Trick . . . . . . . . . . . . . . . . . 7 5.4 O(N log N ) Tandem Repeats . . . . . . . . . . . . . . . . . 20 ll f(ll x){ return a * x + b; }
2.13 O(N ) Smallest Enclosing Circle . . . . . . . . . . . . . . . . 7 5.5 Suffix Automaton . . . . . . . . . . . . . . . . . . . . . . . . 21 };
2.14 O(N + Q log N ) K-D Tree . . . . . . . . . . . . . . . . . . . 7 5.6 Bitset LCS . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 vector<Line> v; int pv;
2.15 O(N log N ) Voronoi Diagram . . . . . . . . . . . . . . . . . 7 5.7 Lyndon Factorization, Minimum Rotation . . . . . . . . . . 21 void init(){ v.clear(); pv = 0; }
5.8 All LCS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 int chk(const Line &a, const Line &b, const Line &c) const {
3 Graph 8 return (__int128_t)(a.b - b.b) * (b.a - c.a) <=
3.1 Euler Tour . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 6 Misc 21 (__int128_t)(c.b - b.b) * (b.a - a.a);
3.2 2-SAT . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 6.1 CMakeLists.txt . . . . . . . . . . . . . . . . . . . . . . . . . 21 }
3.3 Horn SAT . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 6.2 Calendar . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 void insert(Line l){
3.4 2-QBF . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 6.3 Ternary Search . . . . . . . . . . . . . . . . . . . . . . . . . 21 if(v.size() > pv && v.back().a == l.a){ // fix if min query
3.5 BCC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 6.4 Add/Mul Update, Range Sum Query . . . . . . . . . . . . . 21 if(l.b < v.back().b) l = v.back(); v.pop_back();
3.6 Prufer Sequence . . . . . . . . . . . . . . . . . . . . . . . . . 9 6.5 O(N × max W ) Subset Sum (Fast Knapsack) . . . . . . . . 21 }
3.7 Tree Compress . . . . . . . . . . . . . . . . . . . . . . . . . 9 6.6 Monotone Queue Optimization . . . . . . . . . . . . . . . . 21 while(v.size() >= pv+2 && chk(v[v.size()-2], v.back(), l))
3.8 Tree Binarization . . . . . . . . . . . . . . . . . . . . . . . . 9 6.7 Slope Trick . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 v.pop_back();
3.9 O(3V /3 ) Maximal Clique . . . . . . . . . . . . . . . . . . . . 9 6.8 Aliens Trick . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 v.push_back(l);
3.10 O(V log V ) Tree Isomorphism . . . . . . . . . . . . . . . . . 9 6.9 SWAMK, Min Plus Convolution . . . . . . . . . . . . . . . 22 }
3.11 O(E √
log E) Complement Spanning Forest . . . . . . . . . . 9 6.10 Money for Nothing (WF17D) . . . . . . . . . . . . . . . . . 22 p query(ll x){ // if min query, then v[pv].f(x) >= v[pv+1].f(x)
3.12 O(E √ V ) Bipartite Matching, Konig, Dilworth . . . . . . . 10 6.11 Exchange Argument on Tree (WF16L,CERC13E) . . . . . . 22 while(pv+1 < v.size() && v[pv].f(x) <= v[pv+1].f(x)) pv++;
3.13 O(V 2 E) Push Relabel . . . . . . . . . . . . . . . . . . . . 10 6.12 Hook Length Formula . . . . . . . . . . . . . . . . . . . . . 23 return {v[pv].f(x), v[pv].c};
3.14 LR Flow, Circulation . . . . . . . . . . . . . . . . . . . . . . 10 6.13 Floating Point Add (Kahan) . . . . . . . . . . . . . . . . . . 23 }
3.15 Min Cost Circulation . . . . . . . . . . . . . . . . . . . . . . 11 6.14 Random, PBDS, Bit Trick, Bitset . . . . . . . . . . . . . . . 23 ///// line container start (max query) /////
3.16 O(V 3 ) Hungarian Method . . . . . . . . . . . . . . . . . . . 11 6.15 Fast I/O, Fast Div, Fast Mod . . . . . . . . . . . . . . . . . 23 struct Line {
3.17 O(V 3 ) Global Min Cut . . . . . . . . . . . . . . . . . . . . . 11 6.16 DP Optimization . . . . . . . . . . . . . . . . . . . . . . . . 23 mutable ll k, m, p;
3.18 O(V 2 + V√× Flow) Gomory-Hu Tree . . . . . . . . . . . . . 11 6.17 Highly Composite Numbers, Large Prime . . . . . . . . . . 23 bool operator<(const Line& o) const { return k < o.k; }
3.19 O(V + E V ) Count/Find 3/4 Cycle . . . . . . . . . . . . . 11 6.18 DLAS(Diversified Late Acceptance Search) . . . . . . . . . 23 bool operator<(ll x) const { return p < x; }
3.20 O(V log V ) Rectlinear MST . . . . . . . . . . . . . . . . . . 12 6.19 Stack Hack . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 }; // (for doubles, use inf = 1/.0, div(a,b) = a/b)
3.21 O(V E) Shortest Mean Cycle . . . . . . . . . . . . . . . . . 12 6.20 Python Decimal . . . . . . . . . . . . . . . . . . . . . . . . . 24 struct LineContainer : multiset<Line, less<>> {
3.22 O(V 2 ) Stable Marriage Problem . . . . . . . . . . . . . . . 12 6.21 Java I/O . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 static const ll inf = LLONG_MAX; // div: floor
Soongsil University – PS akgwi Page 2 of 25
ll div(ll a, ll b) { return a / b - ((a ^ b) < 0 && a % b); } 1.4 Kinetic Segment Tree void update(int x, const line_t &v, int node=1, int s=0, int
bool isect(iterator x, iterator y) { e=SZ-1){
if (y == end()) return x->p = inf, 0; // 일반적으로 heaten 함수는 교점 s개일 때 O(lambda_{s+2}(n)log^2n) push(node, s, e); int m = (s + e) / 2;
if (x->k == y->k) x->p = x->m > y->m ? inf : -inf; // update가 insert/delete만 주어지면 O(lambda_{s+1}(n)log^2n) if(s == e){ T[node] = v; return; }
else x->p = div(y->m - x->m, x->k - y->k); // update가 없으면 O(lambda_{s}(n)log^2n) if(x <= m)update(x,v, node<<1, s, m), push(node<<1|1, m+1,
return x->p >= y->p; // s = 0: 1 | s = 1: n | s = 2: 2n-1 | s = 3: 2n alpha(n) + O(n) e);
} // s = 4: O(n * 2^alpha(n)) | s = 5: O(n alpha(n) * 2^alpha(n)) else update(x, v, node<<1|1, m+1, e), push(node<<1, s, m);
void add(ll k, ll m) { // apply_heat(heat): x좌표가 heat 증가했을 때의 증가량을 v에 더함 pull(node, s, e);
auto z = insert({k, m, 0}), y = z++, x = y; // heaten(l, r, t): 구간의 온도를 t 만큼 증가 }
while (isect(y, z)) z = erase(z); struct line_t{ void add(int l, int r, ll v, int node=1, int s=0, int e=SZ-1){
if (x != begin() && isect(--x, y)) isect(x, y = erase(y)); ll a, b, v, idx; line_t() : line_t(0, nINF) {} push(node, s, e); int m = (s + e) / 2;
while ((y = x) != begin() && (--x)->p >= y->p) isect(x, line_t(ll a, ll b) : line_t(a, b, -1) {} if(r < s || e < l) return;
erase(y)); line_t(ll a, ll b, ll idx) : a(a), b(b), v(b), idx(idx) {} if(l <= s && e <= r){ T[node].lz_add += v; push(node, s, e);
} ll query(ll x) { assert(!empty()); void apply_heat(ll heat){ v += a * heat; } return; }
auto l = *lower_bound(x); return l.k * x + l.m; } void apply_add(ll lz_add){ v += lz_add; } add(l,r,v,node*2,s,m); add(l,r,v,node*2+1,m+1,e);
}; ll cross(const line_t &l) const { pull(node, s, e);
if(a == l.a) return pINF; ll p = v - l.v, q = l.a - a; }
1.3 Color Processor if(q < 0) p = -p, q = -q; void heaten(int l,int r,ll t,int node=1,int s=0,int e=SZ-1){
template<class CT, class T> struct color_processor { return p >= 0 ? (p + q - 1) / q : -p / q * -1; push(node, s, e); int m = (s + e) / 2;
map<array<CT, 2>, T> v; // CT: coord type } ll cross_after(const line_t &l, ll temp) const { if(r < s || e < l) return;
color_processor(T col={}): v({{{MIN,MAX},col}}){} ll res = cross(l); return res > temp ? res : pINF; } if(l <= s && e <= r){ _heat(t, node, s, e); return; }
auto get_range(CT p){ return *prev(v.upper_bound({p, MAX})); } }; heaten(l,r,t,node*2,s,m); heaten(l,r,t,node*2+1,m+1,e);
// Cover [l, r) with color c, amortized O(1) process call struct range_kinetic_segment_tree{ pull(node, s, e);
// process(l, r, pc, c): color of [l, r) change pc -> c struct node_t{ }
auto cover(CT l, CT r, T c, auto process){ line_t v; ll melt, heat, lz_add; node_t():node_t(line_t()){} void _heat(ll t, int node=1, int s=0, int e=SZ-1){
array<CT, 2> I{l, l}; node_t(ll a, ll b, ll idx) : node_t(line_t(a, b, idx)) {} push(node, s, e); int m = (s + e) / 2;
auto it = v.lower_bound(I); node_t(const line_t &v):v(v),melt(pINF),heat(0),lz_add(0){} if(T[node].melt > t){ T[node].heat += t; push(node, s, e);
if(it != v.begin() && l < prev(it)->fi[1]){ bool operator < (const node_t &o) const { return return; }
auto x = *--it; v.erase(it); tie(v.v,v.a) < tie(o.v.v,o.v.a); } _heat(t,node*2,s,m);_heat(t,node*2+1,m+1,e);pull(node,s,e);
v.insert({{x.fi[0],l}, x.se}); ll cross_after(const node_t &o, ll temp) const { return }
it = v.insert({{l,x.fi[1]}, x.se}).fi; v.cross_after(o.v, temp); } ll query(int l, int r, int node=1, int s=0, int e=SZ-1){
} void apply_lazy(){ v.apply_heat(heat); v.apply_add(lz_add); push(node, s, e); if(r < s || e < l) return nINF;
while(it != v.end() && it->fi[0] < r){ melt -= heat; } if(l <= s && e <= r) return T[node].v.v; int m = (s + e)/2;
if(r < it->fi[1]){ void clear_lazy(){ heat = lz_add = 0; } return max(query(l,r,node<<1,s,m), query(l,r,node<<1|1,m+1,e));
auto x = *it; v.erase(it); void prop_lazy(const node_t &p){ heat += p.heat; lz_add += } // query end
it = v.insert({{x.fi[0],r}, x.se}).fi; p.lz_add; } };
v.insert({{r,x.fi[1]}, x.se}); bool have_lazy() const { return heat != 0 || lz_add != 0; }
} };
process(max(l,it->fi[0]), min(r,it->fi[1]), it->se, c); node_t T[SZ<<1]; range_kinetic_segment_tree(){ clear(); }
I = {min(I[0],it->fi[0]), max(I[1],it->fi[1])}; void clear(){ fill(T, T+SZ*2, node_t()); } 1.5 Lazy LiChao Tree
it = v.erase(it); void pull(int node, int s, int e){
} return v.insert({I, c}); if(s == e) return; /* get_point(x) : get min(f(x)), O(log X)
} const node_t &l = T[node<<1], &r = T[node<<1|1]; range_min(l,r) get min(f(x)), l<=x<=r, O(log X)
// new_color(l, r, pc): return new color for assert(!l.have_lazy() && !r.have_lazy() && insert(l,r,a,b) : insert f(x)=ax+b, l<=x<=r, O(log^2 X)
// [l, r) previous color pc O(Number of color ranges affected) !T[node].have_lazy()); add(l,r,a,b) : add f(x)=ax+b, l<=x<=r, O(log^2 X)
void recolor(CT l, CT r, auto new_color){ T[node] = max(l, r); WARNING: a != 0인 add가 없을 때만 range_min 가능 */
auto left = v.lower_bound({l, l}); T[node].melt = min({ l.melt, r.melt, l.cross_after(r, 0) }); template<typename T, T LE, T RI, T INF=(long long)(4e18)>
if(l < left->fi[0]){ } struct LiChao{
auto [range, c] = *--left; left = v.erase(left); void push(int node, int s, int e){ struct Node{
left = v.insert(left, {{range[0],l},c}); if(!T[node].have_lazy()) return; T[node].apply_lazy(); int l, r; T a, b, mn, aa, bb;
left = v.insert(left, {{l,range[1]},c}); if(s != e) for(auto c : {node<<1, node<<1|1}) Node(){ l = r = 0; a = 0; b = mn = INF; aa = bb = 0; }
} T[c].prop_lazy(T[node]); void apply(){ mn += bb; a += aa; b += bb; aa = bb = 0; }
auto right = v.lower_bound({r, r}); T[node].clear_lazy(); void add_lazy(T _aa, T _bb){ aa += _aa; bb += _bb; }
if(r < right->fi[0]){ } T f(T x) const { return a * x + b; }
auto [range, c] = *--right; right = v.erase(right); void build(const vector<line_t> &lines, int node=1, int s=0, }; vector<Node> seg; LiChao() : seg(2) {}
right = v.insert(right,{{range[0],r},c}); int e=SZ-1){ void make_child(int n){
right = v.insert(right,{{r,range[1]},c}); if(s == e){ T[node] = s < lines.size() ? node_t(lines[s]) : if(!seg[n].l) seg[n].l = seg.size(), seg.emplace_back();
} node_t(); return; } if(!seg[n].r) seg[n].r = seg.size(), seg.emplace_back();
for(auto it=left; it!=right; ++it) int m = (s + e) / 2; }
it->se = new_color(it->fi[0], it->fi[1], it->se); build(lines,node*2,s,m); build(lines,node*2+1,m+1,e); void push(int node, T s, T e){
} pull(node, s, e); if(seg[node].aa || seg[node].bb){
}; } if(s != e){
Soongsil University – PS akgwi Page 3 of 25
make_child(node); void Rotate(){ Splay(x); x->r = nullptr; x->Update();
seg[seg[node].l].add_lazy(seg[node].aa, seg[node].bb); p->Push(); Push(); for(auto y=x; x->p; Splay(x)) y = x->p, Splay(y), y->r = x,
seg[seg[node].r].add_lazy(seg[node].aa, seg[node].bb); if(IsLeft()) r && (r->p = p), p->l = r, r = p; y->Update();
} seg[node].apply(); else l && (l->p = p), p->r = l, l = p; }
} if(!p->IsRoot()) (p->IsLeft() ? p->p->l : p->p->r) = this; int GetDepth(Node *x){Access(x);x->Push();return GetSize(x->l);}
} auto t = p; p = t->p; t->p = this; t->Update(); Update(); Node* GetRoot(Node *x){
void insert(T l, T r, T a, T b, int node=1, T s=LE, T e=RI){ } Access(x);for(x->Push();x->l;x->Push())x=x->l;return Splay(x);
if(r < s || e < l || l > r) return; void Update(){ }
make_child(node); push(node, s, e); T m = (s + e) >> 1; sz = 1 + GetSize(l) + GetSize(r); sum = now + GetSum(l) + Node* GetPar(Node *x){
seg[node].mn=min({seg[node].mn, a*max(s,l)+b,a*min(e,r)+b}); GetSum(r); Access(x); x->Push(); if(!x->l) return nullptr;
if(s < l || r < e){ } x = x->l; for(x->Push(); x->r; x->Push()) x = x->r;
if(l <= m) insert(l, r, a, b, seg[node].l, s, m); void Update(const T &val){ now = val; Update(); } return Splay(x);
if(m+1 <= r) insert(l, r, a, b, seg[node].r, m+1, e); void Push(){ }
return; Update(now + lz); if(flip) swap(l, r); void Link(Node *p, Node *c){ Access(c); Access(p); c->l = p;
} for(auto c : {l, r}) if(c) c->flip ^= flip, c->lz += lz; p->p = c; c->Update(); }
T &sa = seg[node].a, &sb = seg[node].b; lz = 0; flip = false; void Cut(Node *c){ Access(c); c->l->p = nullptr; c->l = nullptr;
if(a*s+b < sa*s+sb) swap(a, sa), swap(b, sb); } c->Update(); }
if(a*e+b >= sa*e+sb) return; }; Node* GetLCA(Node *x, Node *y){
if(a*m+b < sa*m+sb){ Node* rt; Access(x); Access(y); Splay(x); return x->p ? x->p : x;
swap(a,sa); swap(b,sb); insert(l,r,a,b,seg[node].l,s,m); Node* Splay(Node *x, Node *g=nullptr){ }
} else insert(l, r, a, b, seg[node].r, m+1, e); for(g || (rt=x); x->p!=g; x->Rotate()){ Node* Ancestor(Node *x, int k){
} if(!x->p->IsRoot()) x->p->p->Push(); k = GetDepth(x) - k; assert(k >= 0);
void add(T l, T r, T a, T b, int node=1, T s=LE, T e=RI){ x->p->Push(); x->Push(); for(;;x->Push()){
if(r < s || e < l || l > r) return; if(x->p->p != g) (x->IsLeft() ^ x->p->IsLeft() ? x : int s = GetSize(x->l); if(s == k) return Access(x), x;
make_child(node); push(node, s, e); T m = (s + e) >> 1; x->p)->Rotate(); if(s < k) k -= s + 1, x = x->r; else x = x->l;
if(s < l || r < e){ } }
insert(s, m, seg[node].a, seg[node].b, seg[node].l, s, m); x->Push(); return x; }
insert(m+1,e,seg[node].a,seg[node].b,seg[node].r,m+1,e); } void MakeRoot(Node *x){ Access(x); Splay(x); x->flip ^= 1;
seg[node].a = 0; seg[node].b = seg[node].mn = INF; Node* Kth(int k){ x->Push(); }
if(l <= m) add(l, r, a, b, seg[node].l, s, m); for(auto x=rt; ; x=x->r){ bool IsConnect(Node *x, Node *y){return GetRoot(x)==GetRoot(y);}
if(m+1 <= r) add(l, r, a, b, seg[node].r, m+1, e); for(; x->Push(), x->l && x->l->sz > k; x=x->l); void PathUpdate(Node *x, Node *y, T val){
seg[node].mn=min(seg[seg[node].l].mn, if(x->l) k -= x->l->sz; Node *root = GetRoot(x); // original root
seg[seg[node].r].mn); if(!k--) return Splay(x); MakeRoot(x); Access(y); // make x to root, tie with y
return; } Splay(x); x->lz += val; x->Push();
} seg[node].add_lazy(a, b); push(node, s, e); } MakeRoot(root); // Revert
} Node* Gather(int s, int e){ auto t = Kth(e+1); return Splay(t, // edge update without edge vertex...
T get_point(T x, int node=1, T s=LE, T e=RI){ Kth(s-1))->l; } Node *lca = GetLCA(x, y);
if(node == 0) return INF; push(node, s, e); Node* Flip(int s, int e){ auto x = Gather(s, e); x->flip ^= 1; Access(lca); Splay(lca); lca->Push();
T m = (s + e) >> 1, res = seg[node].f(x); return x; } lca->Update(lca->now - val);
if(x <= m) return min(res, get_point(x, seg[node].l, s, m)); Node* Shift(int s, int e, int k){ }
else return min(res, get_point(x, seg[node].r, m+1, e)); if(k >= 0){ // shift to right T VertexQuery(Node *x, Node *y){
} k %= e-s+1; if(k) Flip(s, e), Flip(s, s+k-1), Flip(s+k, e); Node *l = GetLCA(x, y); T ret = l->now;
T range_min(T l, T r, int node=1, T s=LE, T e=RI){ } Access(x); Splay(l); if(l->r) ret = ret + l->r->sum;
if(node == 0 || r < s || e < l || l > r) return INF; else{ // shift to left Access(y); Splay(l); if(l->r) ret = ret + l->r->sum;
push(node, s, e); T m = (s + e) >> 1; k = -k; k %= e-s+1; if(k) Flip(s, e), Flip(s, e-k), return ret;
if(l <= s && e <= r) return seg[node].mn; Flip(e-k+1, e); }
return min({ seg[node].f(max(s,l)), seg[node].f(min(e,r)), } Node* GetQueryResultNode(Node *u, Node *v){
range_min(l, r, seg[node].l, s, m), return Gather(s, e); if(!IsConnect(u, v)) return 0;
range_min(l, r, seg[node].r, m+1, e) }); } MakeRoot(u); Access(v); auto ret = v->l;
} int Idx(Node *x){ return x->l->sz; } while(ret->mx != ret->now){
}; //////////// Link Cut Tree Start //////////// if (ret->l && ret->mx == ret->l->mx) ret = ret->l;
Node* Splay(Node *x){ else ret = ret->r;
1.6 Splay Tree, Link-Cut Tree for(; !x->IsRoot(); x->Rotate()){ }
struct Node{ if(!x->p->IsRoot()) x->p->p->Push(); Access(ret); return ret;
Node *l, *r, *p; bool flip; int sz; T now, sum, lz; x->p->Push(); x->Push(); }
Node(){ l = r = p = nullptr; sz = 1; flip = false; now = sum = if(!x->p->IsRoot()) (x->IsLeft() ^ x->p->IsLeft() ? x :
lz = 0; } x->p)->Rotate();
bool IsLeft() const { return p && this == p->l; } } 2 Geometry
bool IsRoot() const { return !p || (this != p->l && this != x->Push(); return x;
p->r); } 2.1 O(log N ) Point in Convex Polygon
}
friend int GetSize(const Node *x){ return x ? x->sz : 0; } void Access(Node *x){ bool Check(const vector<Point> &v, const Point &pt){
friend T GetSum(const Node *x){ return x ? x->sum : 0; } if(CCW(v[0], v[1], pt) < 0) return false;
Soongsil University – PS akgwi Page 4 of 25
int l = 1, r = v.size() - 1; else{ // down //if(xq < 0) xp=-xp, xq=-xq; if(yq < 0) yp=-yp, yq=-yq //gcd?
while(l < r){ if(sign(CCW(pt, v[m], v[m+1])) < 0) s = m; xp /= xq; yp /= yq;
int m = l + r + 1 >> 1; else if(sign(CCW(pt, v[m], v[s])) < 0) s = m; else e = m; if(s1.x == xp && s1.y == yp) flag |= 1;
if(CCW(v[0], v[m], pt) >= 0) l = m; else r = m - 1; } if(e1.x == xp && e1.y == yp) flag |= 1;
} } if(s2.x == xp && s2.y == yp) flag |= 2;
if(l == v.size() - 1) return CCW(v[0], v.back(), pt) == 0 && if(s && local(pt, v[s-1], v[s], v[s+1])) return s; if(e2.x == xp && e2.y == yp) flag |= 2;
v[0] <= pt && pt <= v.back(); if(e != n && local(pt, v[e-1], v[e], v[e+1])) return e; return {flag ? flag : 4, xp, 1, yp, 1};
return CCW(v[0], v[l], pt) >= 0 && CCW(v[l], v[l+1], pt) >= 0 return -1; }
&& CCW(v[l+1], v[0], pt) >= 0; } P perp() const { return P(-y, x); }
} int Closest(const vector<Point> &v, const Point &out, int now){ #define arg(p, q) atan2(p.cross(q), p.dot(q))
int prv = now > 0 ? now-1 : v.size()-1, nxt = now+1 < v.size() bool circleIntersect(P a,P b,double r1,double r2,pair<P, P>*
2.2 Segment Distance, Segment Reflect ? now+1 : 0, res = now; out){
double Proj(Point a, Point b, Point c){ if(CCW(out, v[now], v[prv]) == 0 && Dist(out, v[res]) > if (a == b) { assert(r1 != r2); return false; }
ll t1 = (b - a) * (c - a), t2 = (a - b) * (c - b); Dist(out, v[prv])) res = prv; P vec = b-a; double d2 = vec.dist2(), sum = r1+r2, dif =
if(t1 * t2 >= 0 && CCW(a, b, c) != 0) if(CCW(out, v[now], v[nxt]) == 0 && Dist(out, v[res]) > r1-r2;
return abs(CCW(a, b, c)) / sqrt(Dist(a, b)); Dist(out, v[nxt])) res = nxt; double p = (d2 + r1*r1 - r2*r2)/(d2*2), h2 = r1*r1 - p*p*d2;
else return 1e18; // INF return res; // if parallel, return closest point to out if (sum*sum < d2 || dif*dif > d2) return false; // use EPS
} } // int point_idx = Closest(convex_hull, pt, plz...
double SegmentDist(Point a[2], Point b[2]){ ConvexTangent(hull + hull[0], pt, +-1) % N); P mid = a + vec*p, per = vec.perp() * sqrt(fmax(0, h2) / d2);
double res = 1e18; // NOTE: need to check intersect ///////// *out = {mid + per, mid - per}; return true;
for(int i=0; i<4; i++) res=min(res,sqrt(Dist(a[i/2],b[i%2]))); double polar(pdd x){ return atan2(x.second, x.first); } }
for(int i=0; i<2; i++) res = min(res, Proj(a[0], a[1], b[i])); int tangent(circle &A, circle &B, pdd des[4]){ // return angle vector<P> circleLine(P c, double r, P a, P b) {
for(int i=0; i<2; i++) res = min(res, Proj(b[0], b[1], a[i])); int top = 0; // outer P ab = b - a, p = a + ab * (c-a) * ab / D2(ab);
return res; double d = size(A.O - B.O), a = polar(B.O - A.O), b = PI + a; T s = (b - a) / (c - a), h2 = r*r - s*s / D2(ab);
} double t = sq(d) - sq(A.r - B.r); if (abs(h2) < EPS) return {p}; if (h2 < 0) return {};
P Reflect(P p1, P p2, P p3){ // line p1-p2, point p3 if (t >= 0){ P h = ab / D1(ab) * sqrtl(h2); return {p - h, p + h};
auto [x1,y1] = p1; auto [x2,y2] = p2; auto [x3,y3] = p3; t = sqrt(t); double p = atan2(B.r - A.r, t); } // use circleLine if you use double...
auto a = y1-y2, b = x2-x1, c = x1 * (y2-y1) + y1 * (x1-x2); des[top++] = pdd(a + p + PI / 2, b + p - PI / 2); int CircleLineIntersect(P o, T r, P p1, P p2, bool segment){
auto d = a * y3 - b * x3; des[top++] = pdd(a - p - PI / 2, b - p + PI / 2); P s = p1, d = p2 - p1; // line : s + kd, int support
T x = -(a*c+b*d) / (a*a+b*b), y = (a*d-b*c) / (a*a+b*b); } T a = d * d, b = (s - o) * d * 2, c = D2(s, o) - r * r;
return 2 * P(x,y) - p3; t = sq(d) - sq(A.r + B.r); // inner T det = b * b - 4 * a * c; // solve ak^2 + bk + c = 0, a > 0
} if (t >= 0){ t = sqrt(t); if(!segment) return Sign(det) + 1;
double p = atan2(B.r + A.r, t); if(det <= 0) return det ? 0 : 0 <= -b && -b <= a + a;
des[top++] = pdd(a + p - PI / 2, b + p - PI / 2); bool f11 = b <= 0 || b * b <= det;
2.3 Tangent Series des[top++] = pdd(a - p + PI / 2, b - p + PI / 2); bool f21 = b <= 0 && b * b >= det;
template<bool UPPER=true> // O(log N) } bool f12 = a+a+b >= 0 && det <= (a+a+b) * (a+a+b);
Point GetPoint(const vector<Point> &hull, real_t slope){ return top; bool f22 = a+a+b >= 0 || det >= (a+a+b) * (a+a+b);
auto chk = [slope](real_t dx, real_t dy){ return UPPER ? dy } return (f11 && f12) + (f21 && f22);
>= slope * dx : dy <= slope * dx; }; pair<T, T> CirclePointTangent(P o, double r, P p){ } // do not use this if you want to use double...
int l = -1, r = hull.size() - 1; T op=D1(p,o), u=atan2l(p.y-o.y, p.x-o.x), v=acosl(r/op); double circlePoly(P c, double r, vector<P> ps){ // return area
while(l + 1 < r){ return {u + v, u - v}; auto tri = [&](P p, P q) { // ps must be ccw polygon
int m = (l + r) / 2; } // COORD 1e4 EPS 1e-7 / COORD 1e3 EPS 1e-9 with circleLine auto r2 = r * r / 2; P d = q - p;
if(chk(hull[m+1].x - hull[m].x, hull[m+1].y - auto a = d.dot(p)/d.dist2(), b = (p.dist2()-r*r)/d.dist2();
hull[m].y)) l = m; else r = m; 2.4 Intersect Series auto det = a * a - b;
} // 0: not intersect, -1: infinity, 4: intersect if (det <= 0) return arg(p, q) * r2;
return hull[r]; // 1/2/3: intersect first/second/both segment corner auto s = max(0., -a-sqrt(det)), t = min(1., -a+sqrt(det));
} // flag, xp, xq, yp, yq : (xp / xq, yp / yq) if (t < 0 || 1 <= s) return arg(p, q) * r2;
int ConvexTangent(const vector<Point> &v, const Point &pt, int using T = __int128_t; // T <= O(COORD^3) P u = p + d * s, v = p + d * t;
up=1){ //given outer point, O(log N) tuple<int,T,T,T,T> SegmentIntersect(P s1, P e1, P s2, P e2){ return arg(p,u) * r2 + u.cross(v)/2 + arg(v,q) * r2;
auto sign = [&](ll c){ return c>0 ? up : c==0 ? 0 : -up; }; if(!Intersect(s1, e1, s2, e2)) return {0, 0, 0, 0, 0}; };
auto local = [&](Point p, Point a, Point b, Point c){ auto det = (e1 - s1) / (e2 - s2); auto sum = 0.0;
return sign(CCW(p, a, b)) <= 0 && sign(CCW(p, b, c)) >= 0; if(!det){ rep(i,0,sz(ps)) sum += tri(ps[i] - c, ps[(i+1)%sz(ps)] - c);
}; // assert(v.size() >= 2); if(s1 > e1) swap(s1, e1); return sum;
int n = v.size() - 1, s = 0, e = n, m; if(s2 > e2) swap(s2, e2); }
if(local(pt, v[1], v[0], v[n-1])) return 0; if(e1 == s2) return {3, e1.x, 1, e1.y, 1}; // extrVertex: point of hull, max projection onto line
while(s + 1 < e){ if(e2 == s1) return {3, e2.x, 1, e2.y, 1}; #define cmp(i,j) sgn(dir.perp().cross(poly[(i)%n]-poly[(j)%n]))
m = (s + e) / 2; return {-1, 0, 0, 0, 0}; #define extr(i) cmp(i + 1, i) >= 0 && cmp(i, i - 1 + n) < 0
if(local(pt, v[m-1], v[m], v[m+1])) return m; } int extrVertex(vector<P>& poly, P dir) {
if(sign(CCW(pt, v[s], v[s+1])) < 0){ // up T p = (s2 - s1) / (e2 - s2), q = det, flag = 0; int n = sz(poly), lo = 0, hi = n;
if(sign(CCW(pt, v[m], v[m+1])) > 0) e = m; T xp = s1.x * q + (e1.x - s1.x) * p, xq = q; if (extr(0)) return 0;
else if(sign(CCW(pt, v[m], v[s])) > 0) s = m; else e = m; T yp = s1.y * q + (e1.y - s1.y) * p, yq = q;
} if(xp%xq || yp%yq) return {4,xp,xq,yp,yq};//gcd?
Soongsil University – PS akgwi Page 5 of 25
while (lo + 1 < hi) { else if(d <= abs(a.r - b.r) + eps) { A += v[j].cross(v[i]);
int m = (lo + hi) / 2; if (extr(m)) return m; if (a.r < b.r) res.emplace_back(0, 2 * pi); } return res / A / 3;
int ls = cmp(lo + 1, lo), ms = cmp(m + 1, m); } else if(d < abs(a.r + b.r) - eps) { }
(ls < ms || (ls == ms && ls == cmp(lo, m)) ? hi : lo) = m; double o = acos((a.r*a.r + d*d - b.r*b.r) / (2 * a.r * d)); // O(points^2), area of union of n polygon, ccw polygon
} double z = NormAngle(atan2((b.o - a.o).y, (b.o - a.o).x)); int sideOf(P s, P e, P p) { return sgn((e-s)/(p-s)); }
return lo; double l = NormAngle(z - o), r = NormAngle(z + o); int sideOf(const P& s, const P& e, const P& p, double eps) {
} if(l > r) res.emplace_back(l, 2*pi), res.emplace_back(0,r); auto a = (e-s)/(p-s); auto l=D1(e-s) * eps;
//(-1,-1): no collision else res.emplace_back(l, r); return (a > l) - (a < -l);
//(i,-1): touch corner } return res; }
//(i,i): along side (i,i+1) }// circle should be identical double rat(P a, P b) { return sgn(b.x) ? a.x/b.x : a.y/b.y; }
//(i,j): cross (i,i+1)and(j,j+1) double CircleUnionArea(vector<Cir> c) { double polyUnion(vector<vector<P>>& poly) { // (points)^2
//(i,i+1): cross corner i int n = c.size(); double a = 0, w; double ret = 0;
// O(log n), ccw no colinear point convex polygon for (int i = 0; w = 0, i < n; ++i) { rep(i,0,sz(poly)) rep(v,0,sz(poly[i])) {
// P perp() const { return P(-y, x); } vector<pair<double, double>> s = {{2 * pi, 9}}, z; P A = poly[i][v], B = poly[i][(v + 1) % sz(poly[i])];
#define cmpL(i) sgn(a.cross(poly[i], b)) for (int j = 0; j < n; ++j) if (i != j) { vector<pair<double, int>> segs = {{0, 0}, {1, 0}};
array<int, 2> lineHull(P a, P b, vector<P>& poly) { // O(log N) z = CoverSegment(c[i], c[j]); rep(j,0,sz(poly)) if (i != j) { // START
int endA = extrVertex(poly, (a - b).perp()); for (auto &e : z) s.push_back(e); } /* for j */ rep(u,0,sz(poly[j])) {
int endB = extrVertex(poly, (b - a).perp()); sort(s.begin(), s.end()); P C = poly[j][u], D = poly[j][(u + 1) % sz(poly[j])];
if (cmpL(endA) < 0 || cmpL(endB) > 0) return {-1, -1}; auto F = [&] (double t) { return c[i].r * (c[i].r * t + int sc = sideOf(A, B, C), sd = sideOf(A, B, D);
array<int, 2> res; c[i].o.x * sin(t) - c[i].o.y * cos(t)); }; if (sc != sd) {
rep(i,0,2) { for (auto &e : s) { double sa = C.cross(D, A), sb = C.cross(D, B);
int lo = endB, hi = endA, n = sz(poly); if (e.first > w) a += F(e.first) - F(w); if (min(sc, sd) < 0)
while ((lo + 1) % n != hi) { w = max(w, e.second); } /* for e */ segs.emplace_back(sa / (sa - sb), sgn(sc - sd));
int m = ((lo + hi + (lo < hi ? 0 : n)) / 2) % n; } return a * 0.5; } }
(cmpL(m) == cmpL(endB) ? lo : hi) = m; else if (!sc && !sd && j<i && sgn((B-A).dot(D-C))>0){
} 2.6 Segment In Polygon segs.emplace_back(rat(C - A, B - A), 1);
res[i] = (lo + !cmpL(hi)) % n; // WARNING: C.push_back(C[0]) before call function segs.emplace_back(rat(D - A, B - A), -1);
swap(endA, endB); bool segment_in_polygon_non_strict(vector<P> &C, P s, P e){ } /*else if*/ } /*rep u*/ } /*rep j*/ // END
} if(!pip(C, s) || !pip(C, e)) return false; sort(all(segs));
if (res[0] == res[1]) return {res[0], -1}; if(s == e) return true; P d = e - s; for (auto& s : segs) s.first = min(max(s.first, 0.0), 1.0);
if (!cmpL(res[0]) && !cmpL(res[1])) vector<pair<frac,int>> v; auto g=raypoints(C, s, d, v); double sum = 0; int cnt = segs[0].second;
switch ((res[0] - res[1] + sz(poly) + 1) % sz(poly)) { for(auto [fr,ev] : v){ // in(06) out(27) rep(j,1,sz(segs)) {
case 0: return {res[0], res[0]}; if(fr.first < 0 || g < fr) continue; if (!cnt) sum += segs[j].first - segs[j - 1].first;
case 2: return {res[1], res[1]}; if(ev == 4) return false; // pass outside corner cnt += segs[j].second;
} if(fr < g && (ev == 2 || ev == 7)) return false; }
return res; if(0 < fr.first && (ev == 0 || ev == 6)) return false; ret += A.cross(B) * sum;
} } return true; } return abs(ret) / 2;
} }
2.5 O(N 2 log N ) Circles Area
ld NormAngle(ld v){while(v < -EPS) v += M_PI * 2; 2.7 Polygon Cut, Center, Union 2.8 Polycon Raycast
while(v > M_PI * 2 + EPS) v -= M_PI * 2; // Returns the polygon on the left of line l // ray A + kd and CCW polygon C, return events {k, event_id}
return v; } // *: dot product, ^: cross product // 0: out->line / 1: in->line / 2: line->out / 3: line->in
ld TwoCircleUnion(const Circle &p, const Circle &q) { // l = p + d*t, l.q() = l + d // 4: pass corner outside / 5: pass corner inside / 6: out -> in
ld d = D1(p.o - q.o); if(d >= p.r+q.r-EPS) return 0; // doubled_signed_area(p,q,r) = (q-p) ^ (r-p) / 7: in -> out
else if(d <= abs(p.r-q.r)+EPS) return pow(min(p.r,q.r),2)*PI; template<class T> vector<point<T>> polygon_cut(const // WARNING: C.push_back(C[0]) before use, ccw, no colinear
ld pc = (p.r*p.r + d*d - q.r*q.r) / (p.r*d*2), pa = acosl(pc); vector<point<T>> &a, const line<T> &l){ struct frac{
ld qc = (q.r*q.r + d*d - p.r*p.r) / (q.r*d*2), qa = acosl(qc); vector<point<T>> res; ll first, second; frac(){}
ld ps = p.r*p.r*pa - p.r*p.r*sin(pa*2)/2; for(auto i = 0; i < (int)a.size(); ++ i){ frac(ll a, ll b) : first(a), second(b) {
ld qs = q.r*q.r*qa - q.r*q.r*sin(qa*2)/2; auto cur = a[i], prev = i ? a[i - 1] : a.back(); if( b < 0 ) first = -a, second = -b; // operator cast int128
return ps + qs; } bool side = doubled_signed_area(l.p, l.q(), cur) > 0; } double v(){ return 1.*first/second; } // operator <,<=,==
ld TwoCircleIntersect(P p1, P p2, ld r1, ld r2){ if(side != (doubled_signed_area(l.p, l.q(), prev) > 0)) }; // assert(d != P(0,0));
auto f = [](ld a, ld b, ld c){ res.push_back(l.p + (cur - l.p ^ prev - cur) / (l.d ^ prev frac raypoints(const vector<P> &C, P A, P d, vector<pair<frac,
return acosl((a*a + b*b - c*c) / (2*a*b)); }; - cur) * l.d); int>> &R){ vector<pair<frac, int>> L;
ld d = D1(p1, p2); if(d + EPS > r1 + r2) return 0; if(side) res.push_back(cur); auto g = gcd(abs(d.x), abs(d.y)); d.x /= g, d.y /= g;
if(d < abs(r1-r2) + EPS) return min(r1,r2)*min(r1,r2)*M_PI; } for(int i = 0; i+1 < C.size(); i++){ P v = C[i+1] - C[i];
ld t1 = f(r1, d, r2), t2 = f(r2, d, r1); return res; int a = sign(d/(C[i]-A)), b = sign(d/(C[i+1]-A));
return r1*r1*(t1-sinl(t1)*cosl(t1)) } if(a == 0)L.emplace_back(frac(d*(C[i]-A)/size2(d), 1), b);
+ r2*r2*(t2-sinl(t2)*cosl(t2)); } P polygonCenter(const vector<P>& v){ // center of mass if(b == 0)L.emplace_back(frac(d*(C[i+1]-A)/size2(d), 1),a);
vector<pair<double, double>> CoverSegment(Cir a, Cir b) { P res(0, 0); double A = 0; if(a*b == -1) L.emplace_back(frac((A-C[i])/v, v/d), 6);
double d = abs(a.o - b.o); vector<pair<double, double>> res; for (int i = 0, j = sz(v) - 1; i < sz(v); j = i++) { } sort(L.begin(), L.end());
if(sign(a.r + b.r - d) == 0); /* skip */ res = res + (v[i] + v[j]) * v[j].cross(v[i]); for(int i = 0; i < L.size(); i++){
Soongsil University – PS akgwi Page 6 of 25
// assert(i+2 >= L.size() || !(L[i].first == L[i+2].first)); return tie(x,f,y) < tie(e.x,e.f,e.y); if(dq.size() && same(CCW(o, dq.back().slope(), i.slope()),
if(i+1<L.size()&&L[i].first==L[i+1].first&&L[i].second!=6){ // strict 0)) continue;
int a = L[i].second, b = L[i+1].second; // return make_tuple(x,-f,y) < make_tuple(e.x,-e.f,e.y); while(dq.size() >= 2 && CheckHPI(dq[dq.size()-2], dq.back(),
R.emplace_back(L[i++].first, a*b? a*b > 0? 4:6:(1-a-b)/2); } i)) dq.pop_back();
} /* end if */ else R.push_back(L[i]); } /* end for */ }; while(dq.size()>=2&&CheckHPI(i,dq[0],dq[1]))dq.pop_front();
int state = 0; // 0: out, 1: in, 2: line+ccw, 3: line+cw tuple<bool,int,int> ShamosHoey(vector<array<Point,2>> v){ dq.push_back(i);
for(auto &[_,n] : R){ int n = v.size(); vector<int> use(n+1); }
if( n == 6 ) n ^= state, state ^= 1; vector<Line> lines; vector<Event> E; multiset<Line> T; while(dq.size() > 2 && CheckHPI(dq[dq.size()-2], dq.back(),
else if( n == 4 ) n ^= state; for(int i=0; i<n; i++){ dq[0])) dq.pop_back();
else if( n == 0 ) n = state, state ^= 2; lines.emplace_back(v[i][0], v[i][1], i); while(dq.size() > 2 && CheckHPI(dq.back(), dq[0], dq[1]))
else if( n == 1 ) n = state^(state>>1), state ^= 3; if(int t=lines[i].get_k(); 0<=t && t<=n) use[t] = 1; dq.pop_front();
} return frac(g, 1); } for(int i=0; i<dq.size(); i++){
} int k = find(use.begin(), use.end(), 0) - use.begin(); Line now = dq[i], nxt = dq[(i+1)%dq.size()];
bool visible(const vector<P> &C, P A, P B){ for(int i=0; i<n; i++){ lines[i].convert_k(k); if(CCW(o, now.slope(), nxt.slope()) <= eps) return
if( A == B ) return true;//return outside? E.emplace_back(lines[i], i, 0); E.emplace_back(lines[i], i, vector<Point>();
char I[4] = "356", O[4] = "157"; 1); ret.push_back(LineIntersect(now, nxt));
vector<pair<frac, int>> R; vector<frac> E; } sort(E.begin(), E.end()); } //for(auto &[x,y] : ret) x = -x, y = -y;
frac s = frac(0, 1), e = raypoints(C, A, B-A, R); for(auto &e : E){ Line::CUR_X = e.x; return ret;
for(auto [f,n] : R){ if(e.f == 0){ } // MakeLine: left side of ray (x1,y1) -> (x2,y2)
if(*find(O, O+3, n+'0')) E.push_back(f); auto it = T.insert(lines[e.i]); Line MakeLine(T x1, T y1, T x2, T y2){
if(*find(I, I+3, n+'0')) E.push_back(f); if(next(it) != T.end() && Intersect(lines[e.i], T a = y2-y1, b = x1-x2, c = x1*a + y1*b; return {a,b,c}; }
} *next(it))) return {true, e.i, next(it)->id};
for(int j = 0; j < E.size(); j += 2) if( !(e <= E[j] || E[j+1] if(it != T.begin() && Intersect(lines[e.i], *prev(it)))
<= s) ) return false; return {true, e.i, prev(it)->id};
2.11 O(M log M ) Dual Graph
return true; } } constexpr int quadrant_id(const Point p){
else{ constexpr int arr[9] = { 5, 4, 3, 6, -1, 2, 7, 0, 1 };
2.9 O(N log N ) Shamos-Hoey auto it = T.lower_bound(lines[e.i]); return arr[sign(p.x)*3+sign(p.y)+4];
struct Line{ if(it != T.begin() && next(it) != T.end() && }
static ll CUR_X; ll x1, y1, x2, y2, id; Intersect(*prev(it), *next(it))) pair<vector<int>, int> dual_graph(const vector<Point> &points,
Line(Point p1, Point p2, int id) : id(id) { return {true, prev(it)->id, next(it)->id}; const vector<pair<int,int>> &edges){
if(p1 > p2) swap(p1, p2); T.erase(it); int n = points.size(), m = edges.size();
tie(x1,y1) = p1; tie(x2,y2) = p2; } vector<int> uf(2*m); iota(uf.begin(), uf.end(), 0);
} Line() = default; } function<int(int)> find = [&](int v){ return v == uf[v] ? v :
int get_k() const { return y1 != y2 ? (x2-x1)/(y1-y2) : -1; } return {false, -1, -1}; uf[v] = find(uf[v]); };
void convert_k(int k){ // x1,y1,x2,y2 = O(COORD^2), use i128 } function<bool(int,int)> merge = [&](int u, int v){ return
in ccw find(u) != find(v) && (uf[uf[u]]=uf[v], true); };
Line res; vector<vector<pair<int,int>>> g(n);
res.x1 = x1 + y1 * k; res.y1 = -x1 * k + y1; 2.10 O(N log N ) Half Plane Intersection for(int i=0; i<m; i++){
res.x2 = x2 + y2 * k; res.y2 = -x2 * k + y2; double CCW(p1, p2, p3); bool same(double a, double b); const g[edges[i].first].emplace_back(edges[i].second, i);
x1 = res.x1; y1 = res.y1; x2 = res.x2; y2 = res.y2; Point o = Point(0, 0); g[edges[i].second].emplace_back(edges[i].first, i);
if(x1 > x2) swap(x1, x2), swap(y1, y2); struct Line{ // ax+by leq c }
} double a, b, c; Line() : Line(0, 0, 0) {} for(int i=0; i<n; i++){
ld get_y(ll offset=0) const { // OVERFLOW Line(double a, double b, double c) : a(a), b(b), c(c) {} const auto base = points[i];
ld t = ld(CUR_X-x1+offset) / (x2-x1); bool operator < (const Line &l) const { sort(g[i].begin(), g[i].end(), [&](auto a, auto b){
return t * (y2 - y1) + y1; } bool f1 = Point(a, b) > o, f2 = Point(l.a, l.b) > o; auto p1=points[a.first]-base, p2=points[b.first]-base;
bool operator < (const Line &l) const { if(f1 != f2) return f1 > f2; return quadrant_id(p1) != quadrant_id(p2) ?
return get_y() < l.get_y(); } double cw = CCW(o, Point(a, b), Point(l.a, l.b)); quadrant_id(p1) < quadrant_id(p2) : p1.cross(p2) > 0;
// strict return same(cw, 0) ? c * hypot(l.a, l.b) < l.c * hypot(a, b) });
/* bool operator < (const Line &l) const { : cw > 0; for(int j=0; j<g[i].size(); j++){
auto le = get_y(), ri = l.get_y(); } Point slope() const { return Point(a, b); } int k = j ? j - 1 : g[i].size() - 1;
if(abs(le-ri) > 1e-7) return le < ri; }; int u = g[i][k].second << 1, v = g[i][j].second << 1 | 1;
if(CUR_X==x1 || CUR_X==l.x1) return get_y(1)<l.get_y(1); Point LineIntersect(Line a, Line b){ auto p1=points[g[i][k].first], p2=points[g[i][j].first];
else return get_y(-1) < l.get_y(-1); double det = a.a*b.b - b.a*a.b, x = (a.c*b.b - a.b*b.c) / det, if(p1 < base) u ^= 1; if(p2 < base) v ^= 1;
} */ y = (a.a*b.c - a.c*b.a) / det; return Point(x, y); merge(u, v);
}; ll Line::CUR_X = 0; } }
struct Event{ // f=0 st, f=1 ed bool CheckHPI(Line a, Line b, Line c){ }
ll x, y, i, f; Event() = default; if(CCW(o, a.slope(), b.slope()) <= 0) return 0; vector<int> res(2*m);
Event(Line l, ll i, ll f) : i(i), f(f) { Point v=LineIntersect(a,b); return v.x*c.a+v.y*c.b>=c.c; for(int i=0; i<2*m; i++) res[i] = find(i);
if(f==0) tie(x,y) = tie(l.x1,l.y1); } auto comp=res;compress(comp);for(auto &i:res)i=IDX(comp,i);
else tie(x,y) = tie(l.x2,l.y2); vector<Point> HPI(vector<Line> v){ int mx_idx = max_element(all(points)) - points.begin();
} sort(v.begin(), v.end()); deque<Line> dq; vector<Point> ret; return {res, res[g[mx_idx].back().second << 1 | 1]};
bool operator < (const Event &e) const { for(auto &i : v){ }
Soongsil University – PS akgwi Page 7 of 25
2.12 O(N 2 log N ) Bulldozer Trick if(x1 <= x && x <= x2) return min((y-y1)*(y-y1), 2.15 O(N log N ) Voronoi Diagram
struct Line{ (y2-y)*(y2-y));
ll i, j, dx, dy; // dx >= 0 if(y1 <= y && y <= y2) return min((x-x1)*(x-x1),
(x2-x)*(x2-x)); /*
Line(int i, int j, const Point &pi, const Point &pj)
T t11 = GetDist(pt, {x1,y1}), t12 = GetDist(pt, {x1,y2}); input: order will be changed, sorted by (y,x) order
: i(i), j(j), dx(pj.x-pi.x), dy(pj.y-pi.y) {}
T t21 = GetDist(pt, {x2,y1}), t22 = GetDist(pt, {x2,y2}); vertex: voronoi intersection points, degree 3, may duplicated
bool operator < (const Line &l) const {
return min({t11, t12, t21, t22}); edge: may contain inf line (-1)
return make_tuple(dy*l.dx, i, j) < make_tuple(l.dy*dx, l.i,
} area
l.j); }
}; - (a,b) = i-th element of area
bool operator == (const Line &l) const {
template<bool IsFirst = 1> struct Cmp { - (u,v) = i-th element of edge
return dy * l.dx == l.dy * dx; }
bool operator() (const Node &a, const Node &b) const { - input[a] is located CCW of u->v line
};
return IsFirst ? a.p.x < b.p.x : a.p.y < b.p.y; } - input[b] is located CW of u->v line
void Solve(){ // V.reserve(N*(N-1)/2)
}; - u->v line is a subset of perpendicular bisector of input[a]
sort(A+1, A+N+1); iota(P+1, P+N+1, 1); vector<Line> V;
struct KDTree { // Warning : no duplicate to input[b] segment
for(int i=1; i<=N; i++) for(int j=i+1; j<=N; j++)
constexpr static size_t NAIVE_THRESHOLD = 16; - Straight line {a, b}, {-1, -1} through midpoint of input[a]
V.emplace_back(i, j, A[i], A[j]);
vector<Node> tree; and input[b]
sort(V.begin(), V.end());
KDTree() = default; */
for(int i=0, j=0; i<V.size(); i=j){
explicit KDTree(const vector<P> &v) { const double EPS = 1e-9;
while(j < V.size() && V[i] == V[j]) j++;
for(int i=0; i<v.size(); i++) tree.emplace_back(v[i], i); int dcmp(double x){ return x < -EPS? -1 : x > EPS ? 1 : 0; }
for(int k=i; k<j; k++){
Build(0, v.size()); // sq(x) = x*x, size(p) = hypot(p.x, p.y)
int u = V[k].i, v = V[k].j; // point id, index -> Pos[id]
} // sz2(p) = sq(p.x)+sq(p.y), r90(p) = (-p.y, p.x)
swap(Pos[u], Pos[v]); swap(A[Pos[u]], A[Pos[v]]);
template<bool IsFirst = 1> double sq(double x){ return x*x; }
if(Pos[u] > Pos[v]) swap(u, v);
void Build(int l, int r) { double size(pdd p){ return hypot(p.x, p.y); }
// @TODO
if(r - l <= NAIVE_THRESHOLD) return; double sz2(pdd p){ return sq(p.x) + sq(p.y); }
}
const int m = (l + r) >> 1; pdd r90(pdd p){ return pdd(-p.y, p.x); }
}
nth_element(tree.begin()+l, tree.begin()+m, tree.begin()+r, pdd line_intersect(pdd a, pdd b, pdd u, pdd v){ return u +
}
Cmp<IsFirst>{}); (((a-u)/b) / (v/b))*v; }
for(int i=l; i<r; i++){ pdd get_circumcenter(pdd p0, pdd p1, pdd p2){
2.13 O(N ) Smallest Enclosing Circle return line_intersect(0.5 * (p0+p1), r90(p0-p1), 0.5 *
tree[m].x1 = min(tree[m].x1, tree[i].p.x); tree[m].y1 =
pt getCenter(pt a, pt b){ return pt((a.x+b.x)/2, (a.y+b.y)/2); } (p1+p2), r90(p1-p2)); }
min(tree[m].y1, tree[i].p.y);
pt getCenter(pt a, pt b, pt c){ double pb_int(pdd left, pdd right, double sweepline){
tree[m].x2 = max(tree[m].x2, tree[i].p.x); tree[m].y2 =
pt aa = b - a, bb = c - a; if(dcmp(left.y - right.y) == 0) return (left.x + right.x) /
max(tree[m].y2, tree[i].p.y);
auto c1 = aa*aa * 0.5, c2 = bb*bb * 0.5, d = aa / bb; 2.0;
}
auto x = a.x + (c1 * bb.y - c2 * aa.y) / d; int sign = left.y < right.y ? -1 : 1;
Build<!IsFirst>(l, m); Build<!IsFirst>(m + 1, r);
auto y = a.y + (c2 * aa.x - c1 * bb.x) / d; pdd v = line_intersect(left, right-left, pdd(0, sweepline),
}
return pt(x, y); } pdd(1, 0));
template<bool IsFirst = 1>
Circle solve(vector<pt> v){ double d1 = sz2(0.5 * (left+right) - v), d2 = sz2(0.5 *
void Query(const P &p, int l, int r, Node &res) const {
pt p = {0, 0}; (left-right));
if(r - l <= NAIVE_THRESHOLD){
double r = 0; int n = v.size(); return v.x + sign * sqrt(std::max(0.0, d1 - d2)); }
for(int i=l; i<r; i++) if(p != tree[i].p && res.dist(p) >
for(int i=0; i<n; i++) if(dst(p, v[i]) > r + EPS){ struct Beachline{
tree[i].dist(p)) res = tree[i];
p = v[i]; r = 0; struct node{ node(){}
}
for(int j=0; j<i; j++) if(dst(p, v[j]) > r + EPS){ node(pdd point, int idx):point(point), idx(idx), end(0),
else{ // else 1
p = getCenter(v[i], v[j]); r = dst(p, v[i]); link{0, 0}, par(0), prv(0), nxt(0) {}
const int m = (l + r) >> 1;
for(int k=0; k<j; k++) if(dst(p, v[k]) > r + EPS){ pdd point; int idx; int end;
const T t = IsFirst ? p.x - tree[m].p.x : p.y -
p = getCenter(v[i], v[j], v[k]); r = dst(v[k], p); node *link[2], *par, *prv, *nxt; };
tree[m].p.y;
}}} node *root;
if(p != tree[m].p && res.dist(p) > tree[m].dist(p)) res =
return {p, r}; } double sweepline;
tree[m];
if(!tree[m].contain(p) && tree[m].dist_to_border(p) >= Beachline() : sweepline(-1e20), root(NULL){ }
2.14 O(N + Q log N ) K-D Tree res.dist(p)) return; inline int dir(node *x){ return x->par->link[0] != x; }
T GetDist(const P &a, const P &b){ return (a.x-b.x) * (a.x-b.x) if(t < 0){ void rotate(node *n){
+ (a.y-b.y) * (a.y-b.y); } Query<!IsFirst>(p, l, m, res); node *p = n->par; int d = dir(n);
struct Node{ if(t*t < res.dist(p)) Query<!IsFirst>(p, m+1, r, res); p->link[d] = n->link[!d];
P p; int idx; } if(n->link[!d]) n->link[!d]->par = p;
T x1, y1, x2, y2; else{ // else 2 n->par = p->par; if(p->par) p->par->link[dir(p)] = n;
Node(const P &p, const int idx) : p(p), idx(idx), x1(1e9), Query<!IsFirst>(p, m+1, r, res); n->link[!d] = p; p->par = n;
y1(1e9), x2(-1e9), y2(-1e9) {} if(t*t < res.dist(p)) Query<!IsFirst>(p, l, m, res); } void splay(node *x, node *f = NULL){
bool contain(const P &pt)const{ return x1 <= pt.x && pt.x <= } /*else 1*/ } /*else 2*/ } /*void Query*/ while(x->par != f){
x2 && y1 <= pt.y && pt.y <= y2; } int Query(const P& p) const { if(x->par->par == f);
T dist(const P &pt) const { return idx == -1 ? INF : Node ret(make_pair<T>(1e9, 1e9), -1); Query(p, 0, else if(dir(x) == dir(x->par)) rotate(x->par);
GetDist(p, pt); } tree.size(), ret); return ret.idx; } else rotate(x);
T dist_to_border(const P &pt) const { }; rotate(x); }
const auto [x,y] = pt; if(f == NULL) root = x;
Soongsil University – PS akgwi Page 8 of 25
} void insert(node *n, node *p, int d){ auto write_edge = [&](int idx, int v){ idx%2 == 0 ? void Get(){
splay(p); node* c = p->link[d]; edge[idx/2].x = v : edge[idx/2].y = v; }; for(int i=1; i<=pv; i++) if(In[i] < Out[i]){ DFS(i); return; }
n->link[d] = c; if(c) c->par = n; auto add_event = [&](BNode* cur){ double nxt; for(int i=1; i<=pv; i++) if(Out[i]){ DFS(i); return; }
p->link[d] = n; n->par = p; if(bl.get_event(cur, nxt)) events.emplace(nxt, cur); }; }// WARNING: path.size() == M + 1 && not trail
node *prv = !d?p->prv:p, *nxt = !d?p:p->nxt; int n = input.size(), cnt = 0;
n->prv = prv; if(prv) prv->nxt = n; arr = new BNode[n*4]; sz = 0; 3.2 2-SAT
n->nxt = nxt; if(nxt) nxt->prv = n; sort(input.begin(), input.end(), [](const pdd &l, const pdd
} void erase(node* n){ &r){ int SZ; vector<vector<int>> G1, G2;
node *prv = n->prv, *nxt = n->nxt; return l.y != r.y ? l.y < r.y : l.x < r.x; }); void Init(int n){ SZ = n; G1 = G2 = vector<vector<int>>(SZ*2); }
if(!prv && !nxt){ if(n == root) root = NULL; return; } BNode* tmp = bl.root = new_node(input[0], 0), *t2; int New(){
n->prv = NULL; if(prv) prv->nxt = nxt; for(int i = 1; i < n; i++){ for(int i=0;i<2;i++) G1.emplace_back(), G2.emplace_back();
n->nxt = NULL; if(nxt) nxt->prv = prv; if(dcmp(input[i].y - input[0].y) == 0){ return SZ++; }
splay(n); add_edge(-1, -1, i-1, i, 0, tmp); void AddEdge(int s, int e){ G1[s].push_back(e);
if(!nxt){ bl.insert(t2 = new_node(input[i], i), tmp, 1); G2[e].push_back(s); }
root->par = NULL; n->link[0] = NULL; tmp = t2; // T(x) = x << 1, F(x) = x << 1 | 1, I(x) = x ^ 1
root = prv; } } void AddCNF(int a, int b){ AddEdge(I(a), b); AddEdge(I(b), a); }
else{ else events.emplace(input[i].y, i); void MostOne(vector<int> vec){ compress(vec);
splay(nxt, n); node* c = n->link[0]; } for(int i=0; i<vec.size(); i++){
nxt->link[0] = c; c->par = nxt; n->link[0] = NULL; while(events.size()){ int now = New();
n->link[1] = NULL; nxt->par = NULL; event q = events.top(); events.pop(); AddEdge(vec[i], T(now)); AddEdge(F(now), I(vec[i]));
root = nxt; } BNode *prv, *cur, *nxt, *site; if(i == 0) continue;
} bool get_event(node* cur, double &next_sweep){ int v = vertex.size(), idx = q.idx; AddEdge(T(now-1), T(now)); AddEdge(F(now), F(now-1));
if(!cur->prv || !cur->nxt) return false; bl.sweepline = q.sweep; AddEdge(T(now-1), I(vec[i])); AddEdge(vec[i], F(now-1));
pdd u = r90(cur->point - cur->prv->point); if(q.type == 0){ }
pdd v = r90(cur->nxt->point - cur->point); pdd point = input[idx]; }
if(dcmp(u/v) != 1) return false; cur = bl.find_bl(point.x);
pdd p = get_circumcenter(cur->point, cur->prv->point, bl.insert(site = new_node(point, idx), cur, 0); 3.3 Horn SAT
cur->nxt->point); bl.insert(prv = new_node(cur->point, cur->idx), site, 0); /* n : numer of variance
next_sweep = p.y + size(p - cur->point); return true; add_edge(-1, -1, cur->idx, idx, site, prv); {}, 0 : x1 | {0, 1}, 2 : (x1 and x2) => x3, (-x1 or -x2 or x3)
} node* find_bl(double x){ add_event(prv); add_event(cur); fail -> empty vector */
node* cur = root; } vector<int> HornSAT(int n, const vector<vector<int>> &cond,
while(cur){ else{ const vector<int> &val){
double left = cur->prv ? pb_int(cur->prv->point, cur = q.cur, prv = cur->prv, nxt = cur->nxt; int m = cond.size(); vector<int> res(n), margin(m), stk;
cur->point, sweepline) : -1e30; if(!prv || !nxt || prv->idx != q.prv || nxt->idx != q.nxt) vector<vector<int>> gph(n);
double right = cur->nxt ? pb_int(cur->point, continue; for(int i=0; i<m; i++){
cur->nxt->point, sweepline) : 1e30; vertex.push_back(get_circumcenter(prv->point, nxt->point, margin[i] = cond[i].size();
if(left <= x && x <= right){ splay(cur); return cur; } cur->point)); if(cond[i].empty()) stk.push_back(i);
cur = cur->link[x > right]; } write_edge(prv->end, v); write_edge(cur->end, v); for(auto j : cond[i]) gph[j].push_back(i);
} add_edge(v, -1, prv->idx, nxt->idx, 0, prv); }
}; using BNode = Beachline::node; bl.erase(cur); while(!stk.empty()){
static BNode* arr; add_event(prv); add_event(nxt); int v = stk.back(), h = val[v]; stk.pop_back();
static int sz; } if(h < 0) return vector<int>();
static BNode* new_node(pdd point, int idx){ } if(res[h]) continue; res[h] = 1;
arr[sz] = BNode(point, idx); return arr + (sz++); } delete arr; for(auto i : gph[h]) if(!--margin[i]) stk.push_back(i);
struct event{ } } return res;
event(double sweep, int idx):type(0), sweep(sweep), idx(idx){} }
event(double sweep, BNode* cur):type(1), sweep(sweep), 3 Graph
prv(cur->prv->idx), cur(cur), nxt(cur->nxt->idx){}
int type, idx, prv, nxt; BNode* cur; double sweep;
3.1 Euler Tour 3.4 2-QBF
bool operator>(const event &l)const{ return sweep > l.sweep; } // Not Directed / Cycle // con[i] \in {A(∀), E(∃)}, 0-based string
}; constexpr int SZ = 1010; // variable: 1-based(parameter), 0-based(computing)
void VoronoiDiagram(vector<pdd> &input, vector<pdd> &vertex, int N, G[SZ][SZ], Deg[SZ], Work[SZ]; // (a or not b) -> {a, -b} in 1-based index
vector<pii> &edge, vector<pii> &area){ void DFS(int v){ // return empty vector if satisfiable, else any solution
Beachline bl = Beachline(); for(int &i=Work[v]; i<=N; i++) while(G[v][I]) G[v][i]--, // T(x) = x << 1, F(x) = x << 1 | 1, I(x) = x ^ 1
priority_queue<event, vector<event>, greater<event>> events; G[i][v]--, DFS(i); vector<int> TwoQBF(int n, string con, vector<pair<int,int>>
auto add_edge = [&](int u, int v, int a, int b, BNode* c1, cout << v << " "; cnf){
BNode* c2){ } auto f = [](int v){ return v > 0 ? T(v-1) : F(-v-1); };
if(c1) c1->end = edge.size()*2; // Directed / Path for(auto &[a,b] : cnf) AddCNF(a=f(a), b=f(b));
if(c2) c2->end = edge.size()*2 + 1; void DFS(int v){ if(!TwoSAT(n)) return {}; int k = SCC.size();
edge.emplace_back(u, v); area.emplace_back(a, b); for(int i=1; i<=pv; i++) while(G[v][i]) G[v][i]--, DFS(i); vector<int> has(k,-1), from(k), to(k), res(n, -1);
}; Path.push_back(v); for(int i=n-1; i>=0; i--){ // WARNING: index is scc id
} if(has[C[T(i)]] != -1 || has[C[F(i)]] != -1) return {};
Soongsil University – PS akgwi Page 9 of 25
if(con[i] == 'A') has[C[T(i)]] = T(i), has[C[F(i)]] = F(i); for(int i=1; i<=n; i++) sort(G[i].begin(), G[i].end()); 3.9 O(3V /3 ) Maximal Clique
} for(int i=1; i<=n; i++) if(P[i] == -1) dfs(i); using B = bitset<128>; template<typename F> //0-based
for(int i=0; i<k; i++) if(has[i] != -1) from[i] = to[i] = 1; return move(res); // sort(all(res)); void maximal_cliques(vector<B>&g,F f,B P=~B(),B X={},B R={}){
for(int i=0; i<n+n; i++){ } if(!P.any()){ if(!X.any()) f(R); return; }
for(auto j : Gph[i]) if(C[i] != C[j]) vector<int> BCC[MAX_V]; // BCC[v] = components which contains v auto q = (P|X)._Find_first(); auto c = P & ~g[q];
DAG[C[i]].push_back(C[j]); void TwoConnected(int n){ // allow multi edge, no self loop for(int i=0; i<g.size(); i++) if(c[i]) {
} int cnt = 0; array<char,MAX_V> vis; vis.fill(0); R[i] = 1; cliques(g, f, P&g[i], X&g[i], R);
for(int i=k-1; i>=0; i--){ function<void(int,int)> dfs = [&dfs,&vis,&cnt](int v, int c){ R[i]=P[i]=0; X[i] = 1; } // faster for sparse gph
bool flag = false; // i -> A? vis[v] = 1; if(c > 0) BCC[v].push_back(c); } // undirected, self loop not allowed, O(3^{n/3})
for(auto j : DAG[i]) flag |= to[j]; for(auto i : G[v]){ B max_independent_set(vector<vector<int>> g){ //g=adj matrix
if(flag && to[i]) return {}; to[i] |= flag; if(vis[i]) continue; int n = g.size(), i, j; vector<B> G(n); B res{};
} if(In[v] <= Low[i]) BCC[v].push_back(++cnt), dfs(i, cnt); auto chk_mx = [&](B a){ if(a.count()>res.count()) res=a; };
for(int i=0; i<k; i++) for(auto j : DAG[i]) from[j] |= else dfs(i, c); for(i=0; i<n; i++) for(int j=0; j<n; j++)
from[i]; // A->i? }}; if(i!=j && !g[i][j])G[i][j]=1;
for(int i=0; i<k; i++){ for(int i=1; i<=n; i++) if(!vis[i]) dfs(i, 0); cliques(G, chk_mx); return res; }
if(has[i] != -1) for(auto v : SCC[i]) res[v/2] = v % 2 == for(int i=1; i<=n; i++) if(BCC[i].empty())
has[i] % 2; BCC[i].push_back(++cnt);
else if(from[i]) for(auto v : SCC[i]) res[v/2] = v % 2 == 0; }void TwoEdgeConnected(int n){}//remove bridge, do flood fill 3.10 O(V log V ) Tree Isomorphism
else if(to[i]) for(auto v : SCC[i]) res[v/2] = v % 2 == 1;
struct Tree{ // (M1,M2)=(1e9+7, 1e9+9), P1,P2 = random int
} 3.6 Prufer Sequence array(sz >= N+2)
for(int i=0; i<n; i++) if(res[i]==-1) res[i]=C[F(i)]<C[T(i)];
vector<pair<int,int>> Gen(int n, vector<int> a){ // a : int N; vector<vector<int>> G; vector<pair<int,int>> H;
return res;
[1,n]^(n-2) vector<int> S, C; // size,centroid
}
if(n == 1) return {}; if(n == 2) return { make_pair(1, 2) }; Tree(int N) : N(N), G(N+2), H(N+2), S(N+2) {}
3.5 BCC vector<int> deg(n+1); for(auto i : a) deg[i]++; void addEdge(int s, int e){ G[s].push_back(e);
// Call tarjan(N) before use!!! vector<pair<int,int>> res; priority_queue<int> pq; G[e].push_back(s); }
vector<int> G[MAX_V]; int In[MAX_V], Low[MAX_V], P[MAX_V]; for(int i=n; i; i--) if(!deg[i]) pq.emplace(i); int getCentroid(int v, int b=-1){
void addEdge(int s,int e){G[s].push_back(e);G[e].push_back(s);} for(auto i : a){ S[v] = 1; // do not merge if-statements
void tarjan(int n){ /// Pre-Process res.emplace_back(i, pq.top()); pq.pop(); for(auto i : G[v]) if(i!=b) if(int now=getCentroid(i,v);
int pv = 0; if(!--deg[i]) pq.push(i); now<=N/2) S[v]+=now; else break;
function<void(int,int)> dfs = [&pv,&dfs](int v, int b){ }int u = pq.top(); pq.pop(); int v = pq.top(); pq.pop(); if(N - S[v] <= N/2) C.push_back(v); return S[v] = S[v];
In[v] = Low[v] = ++pv; P[v] = b; res.emplace_back(u, v); return res; }
for(auto i : G[v]){ } int init(){
if(i == b) continue; getCentroid(1); if(C.size() == 1) return C[0];
if(!In[i]) dfs(i, v), Low[v] = min(Low[v], Low[i]); else 3.7 Tree Compress int u = C[0], v = C[1], add = ++N;
Low[v] = min(Low[v], In[i]); G[u].erase(find(G[u].begin(), G[u].end(), v));
void Go(int v, const vector<int> vertex, G[v].erase(find(G[v].begin(), G[v].end(), u));
}}; vector<tuple<int,int,int>> &edge){
for(int i=1; i<=n; i++) if(!In[i]) dfs(i, -1); G[add].push_back(u); G[u].push_back(add);
static int pv = 1; // vertex is sorted by dfs order G[add].push_back(v); G[v].push_back(add);
} while(pv < vertex.size() && In[vertex[pv]] <= Out[v]){
vector<int> cutVertex(int n){ return add;
int nxt = vertex[pv++]; }
vector<int> res; array<char,MAX_V> isCut; isCut.fill(0); edge.emplace_back(v, nxt, D[nxt]-D[v]);
function<void(int)> dfs = [&dfs,&isCut](int v){ pair<int,int> build(const vector<ll> &P1, const vector<ll>
Go(nxt, vertex, edge); &P2, int v, int b=-1){
int ch = 0; }}
for(auto i : G[v]){ vector<pair<int,int>> ch; for(auto i : G[v]) if(i != b)
if(P[i] != v) continue; dfs(i); ch++; ch.push_back(build(P1, P2, i, v));
if(P[v] == -1 && ch > 1) isCut[v] = 1; 3.8 Tree Binarization ll h1 = 0, h2 = 0; stable_sort(ch.begin(), ch.end());
else if(P[v] != -1 && Low[i] >= In[v]) isCut[v]=1; void make_binary(int v=1, int real=1, int b=-1, int idx=0){ if(ch.empty()){ return {1, 1}; }
}}; for(int t=idx; t<Inp[v].size(); t++){ for(int i=0; i<ch.size(); i++)
for(int i=1; i<=n; i++) if(P[i] == -1) dfs(i); auto i = Inp[v][t]; if(i.first == b) continue; h1=(h1+(ch[i].first^P1[P1.size()-1-i])*P1[i])%M1,
for(int i=1; i<=n; i++) if(isCut[i]) res.push_back(i); if(G[real].empty() || t+1 == Inp[v].size() h2=(h2+(ch[i].second^P2[P2.size()-1-i])*P2[i])%M2;
return move(res); || t+2 == Inp[v].size() && Inp[v][t+1].first == b){ return H[v] = {h1, h2};
} G[real].push_back(i);//do not change order!!! }
vector<PII> cutEdge(int n){ make_binary(i.first, i.first, v); int build(const vector<ll> &P1, const vector<ll> &P2){
vector<PII> res; G[i.first].emplace_back(real, i.second); int rt = init(); build(P1, P2, rt); return rt;
function<void(int)> dfs = [&dfs,&res](int v){ continue; }
for(int t=0; t<G[v].size(); t++){ } int nxt = N + ++pv;;//do not change order!!! };
int i = G[v][t]; if(t != 0 && G[v][t-1] == G[v][t]) G[real].emplace_back(nxt, 0);
continue; make_binary(v, nxt, b, t);
if(P[i] != v) continue; dfs(i); G[nxt].emplace_back(real, 0); 3.11 O(E log E) Complement Spanning Forest
if((t+1 == G[v].size() || i != G[v][t+1]) && Low[i] > break; vector<pair<int,int>> ComplementSpanningForest(int n, const
In[v]) res.emplace_back(min(v,i), max(v,i)); } vector<pair<int,int>> &edges){ // V+ElgV
}}; // sort edges if multi edge exist } vector<vector<int>> g(n);
Soongsil University – PS akgwi Page 10 of 25
for(const auto &[u,v] : edges) g[u].push_back(v), } void addEdge(int s, int e, flow_t x){
g[v].push_back(u); tuple<vector<int>, vector<int>, int> minimum_vertex_cover(){ g[s].emplace_back(s, e, x, (int)g[e].size());
for(int i=0; i<n; i++) sort(g[i].begin(), g[i].end()); int matching = maximum_matching(); vector<int> lv, rv; if(s == e) g[s].back().r++;
set<int> alive; fill(track.begin(), track.end(), 0); g[e].emplace_back(e, s, 0, (int)g[s].size()-1);
for(int i=0; i<n; i++) alive.insert(i); for(int i=0; i<n; i++) if(le[i] == -1) dfs_track(i); }
vector<pair<int,int>> res; for(int i=0; i<n; i++) if(!track[i]) lv.push_back(i); void enqueue(int v){
while(!alive.empty()){ for(int i=0; i<m; i++) if(track[n+i]) rv.push_back(i); if(!active[v] && excess[v] > 0 && dist[v] < n){
int u = *alive.begin(); alive.erase(alive.begin()); return {lv, rv, lv.size() + rv.size()}; // s(lv)+s(rv)=mat active[v] = true; bucket[dist[v]].push_back(v); b = max(b,
queue<int> que; que.push(u); } dist[v]); }
while(!que.empty()){ tuple<vector<int>, vector<int>, int> }
int v = que.front(); que.pop(); maximum_independent_set(){ void push(edge_t &e){
for(auto it=alive.begin(); it!=alive.end(); ){ auto [a,b,matching] = minimum_vertex_cover(); flow_t fl = min(excess[e.u], e.c - e.f);
if(auto t=lower_bound(g[v].begin(), g[v].end(), *it); t vector<int> lv, rv; lv.reserve(n-a.size()); if(dist[e.u] == dist[e.v] + 1 && fl > flow_t(0)){
!= g[v].end() && *it == *t) ++it; rv.reserve(m-b.size()); e.f += fl; g[e.v][e.r].f -= fl; excess[e.u] -= fl;
else que.push(*it), res.emplace_back(u, *it), it = for(int i=0, j=0; i<n; i++){ excess[e.v] += fl; enqueue(e.v); }
alive.erase(it); while(j < a.size() && a[j] < i) j++; }
}}}return res; if(j == a.size() || a[j] != i) lv.push_back(i); void gap(int k){
} } for(int i=0; i<n; i++){
for(int i=0, j=0; i<m; i++){ if(dist[i] >= k) count[dist[i]]--, dist[i] = max(dist[i],
√ while(j < b.size() && b[j] < i) j++; n), count[dist[i]]++;
3.12 O(E V ) Bipartite Matching, Konig, Dilworth if(j == b.size() || b[j] != i) rv.push_back(i); enqueue(i); }
struct HopcroftKarp{ } // s(lv)+s(rv)=n+m-mat }
int n, m; vector<vector<int>> g; return {lv, rv, lv.size() + rv.size()}; void relabel(int v){
vector<int> dst, le, ri; vector<char> visit, track; } count[dist[v]]--; dist[v] = n;
HopcroftKarp(int n, int m) : n(n), m(m), g(n), dst(n), le(n, vector<vector<int>> minimum_path_cover(){ // n == m for(const auto &e : g[v]) if(e.c - e.f > 0) dist[v] =
-1), ri(m, -1), visit(n), track(n+m) {} int matching = maximum_matching(); min(dist[v], dist[e.v] + 1);
void add_edge(int s, int e){ g[s].push_back(e); } vector<vector<int>> res; res.reserve(n - matching); count[dist[v]]++; enqueue(v);
bool bfs(){ bool res = false; queue<int> que; fill(track.begin(), track.end(), 0); }
fill(dst.begin(), dst.end(), 0); auto get_path = [&](int v) -> vector<int> { void discharge(int v){
for(int i=0; i<n; i++)if(le[i] == -1)que.push(i),dst[i]=1; vector<int> path{v}; // ri[v] == -1 for(auto &e : g[v]) if(excess[v] > 0) push(e); else break;
while(!que.empty()){ int v = que.front(); que.pop(); while(le[v] != -1) path.push_back(v=le[v]); if(excess[v] > 0) if(count[dist[v]] == 1) gap(dist[v]);
for(auto i : g[v]){ return path; else relabel(v);
if(ri[i] == -1) res = true; }; }
else if(!dst[ri[i]])dst[ri[i]]=dst[v]+1,que.push(ri[i]); for(int i=0; i<n; i++) if(!track[n+i] && ri[i] == -1) flow_t maximumFlow(int _n, int s, int t){
} res.push_back(get_path(i)); // memset dist, excess, count, active 0
} return res; // sz(res) = n-mat n = _n; b = 0; for(auto &e : g[s]) excess[s] += e.c;
return res; } count[s] = n; enqueue(s); active[t] = true;
} vector<int> maximum_anti_chain(){ // n = m while(b >= 0){
bool dfs(int v){ auto [a,b,matching] = minimum_vertex_cover(); if(bucket[b].empty()) b--;
if(visit[v]) return false; visit[v] = 1; vector<int> res; res.reserve(n - a.size() - b.size()); else{
for(auto i : g[v]){ for(int i=0, j=0, k=0; i<n; i++){ int v = bucket[b].back(); bucket[b].pop_back();
if(ri[i] == -1 || !visit[ri[i]] && dst[ri[i]] == dst[v] + while(j < a.size() && a[j] < i) j++; active[v] = false; discharge(v);
1 && dfs(ri[i])){ le[v] = i; ri[i] = v; return true; } while(k < b.size() && b[k] < i) k++; } /*else*/ } /*while*/ return excess[t];
} return false; if((j == a.size() || a[j] != i) && (k == b.size() || b[k] }
} != i)) res.push_back(i); };
int maximum_matching(){ } return res; // sz(res) = n-mat
int res = 0; fill(all(le), -1); fill(all(ri), -1); }
while(bfs()){ };
fill(visit.begin(), visit.end(), 0);
3.14 LR Flow, Circulation
for(int i=0; i<n; i++) if(le[i] == -1) res += dfs(i); 2
√ struct LR_Flow{ LR_Flow(int n) : F(n+2), S(0) {}
} return res; 3.13 O(V E) Push Relabel void add_edge(int s, int e, flow_t l, flow_t r){
} template<typename flow_t> struct Edge { S += abs(l); F.add_edge(s+2, e+2, r-l);
vector<pair<int,int>> maximum_matching_edges(){ int u, v, r; flow_t c, f; Edge() = default; if(l > 0) F.add_edge(s+2, 1, l), F.add_edge(0, e+2, l);
int matching = maximum_matching(); Edge(int u, int v, flow_t c, int r) : u(u), v(v), r(r), c(c), else F.add_edge(0, s+2, -l), F.add_edge(e+2, 1, -l);
vector<pair<int,int>> edges; edges.reserve(matching); f(0) {} } Dinic<flow_t, MAX_U> F; flow_t S;
for(int i=0; i<n; i++) if(le[i] != -1) edges.emplace_back(i, }; bool solve(int s, int t){//maxflow: run F.maximum_flow(s, t)
le[i]); template<typename flow_t, size_t _Sz> struct PushRelabel { if(s != -1) F.add_edge(t+2, s+2, MAX_U); //min cost circ
return edges; using edge_t = Edge<flow_t>; return F.maximum_flow(0,1) == S; }
} int n, b, dist[_Sz], count[_Sz+1]; flow_t get_flow(int s, int e) const { s += 2; e += 2;
void dfs_track(int v){ flow_t excess[_Sz]; bool active[_Sz]; for(auto i : F.g[s]) if(i.c > 0 && i.v == e) return i.f; }
if(track[v]) return; track[v] = 1; vector<edge_t> g[_Sz]; vector<int> bucket[_Sz]; }; struct Circulation{ // demand[i] = in[i] - out[i]
for(auto i : g[v]) track[n+i] = 1, dfs_track(ri[i]); void clear(){ for(int i=0; i<_Sz; i++) g[i].clear(); } Circulation(int n, const vector<flow_t> &demand) : F(n+2) {
Soongsil University – PS akgwi Page 11 of 25
// demand[i] > 0: add_edge(0, i+2, demand[i], demand[i]) if(excess[x] == 0) break; relabel(x); }
// demand[i] <= 0: add_edge(i+2, 1, -demand[i], demand[i]) } /* excess end */ } /* que end */ } return {mn, mn_cut};
} LR_Flow<flow_t, MAX_U> F; } /* while true end */ T ans = 0; }
void add_edge(int s, int e, flow_t l, flow_t r){ for(int i=0; i<n; i++) for(auto &j : g[i])
F.add_edge(s+2, e+2, l, r); } j.cost /= n + 1, ans += j.cost * (j.cap - j.rem);
bool feasible(){ return F.feasible(0, 1); } }; return make_pair(true, ans / 2); 3.18 O(V 2 + V × Flow) Gomory-Hu Tree
} //0-based, S-T cut in graph=S-T cut in gomory-hu tree (path min)
3.15 Min Cost Circulation void bellmanFord(){ vector<Edge> GomoryHuTree(int n, const vector<Edge> &e){
template <class T> struct MinCostCirculation { fill(p.begin(), p.end(), T(0)); bool upd = 1; Dinic<int,100> Flow; vector<Edge> res(n-1); vector<int> pr(n);
const int SCALE = 3; // scale by 1/(1 << SCALE) while(upd){ upd = 0; for(int i=1; i<n; i++, Flow.clear()){ // // bi-directed edge
const T INF = numeric_limits<T>::max() / 2; for(int i = 0; i < n; i++) for(auto &j : g[i]) for(const auto &[s,e,x] : e) Flow.AddEdge(s, e, x);
struct EdgeStack { int s, e; T l, r, cost; }; if(j.rem > 0 && p[j.pos] > p[i] + j.cost) p[j.pos] = int fl = Flow.MaxFlow(pr[i], i);
struct Edge { int pos, rev; T rem, cap, cost; }; p[i] + j.cost, upd = 1; for(int j=i+1; j<n; j++){
int n; vector<vector<Edge>> g; vector<EdgeStack> estk; } if(!Flow.Level[i] == !Flow.Level[j] && pr[i] == pr[j])
LR_Flow<T, 1LL<<60> circ; vector<T> p; } pr[j] = i;
MinCostCirculation(int k) : n(k), g(k), circ(k), p(k) {} }; } /*for-j end*/ res[i-1] = Edge(pr[i], i, fl);
void add_edge(int s, int e, T l, T r, T cost){
} /*for-i end*/ return res; }
estk.push_back({s, e, l, r, cost}); }
3.16 O(V 3 ) Hungarian Method
pair<bool, T> solve(){
for(auto &i:estk)if(i.s!=i.e)circ.add_edge(i.s,i.e,i.l,i.r); // C[j][w] = cost(j-th job, w-th worker), j <= w, O(J^2W) √
if(!circ.solve(-1, -1)) return make_pair(false, T(0)); // ret[i] = minimum cost to assign 0..i jobs to distinct workers 3.19 O(V + E V ) Count/Find 3/4 Cycle
vector<int> cnt(n); T eps = 0; template<typename T>bool ckmin(T &a, const T &b){return b<a ? vector<tuple<int,int,int>> Find3Cycle(int n, const
for(auto &i : estk){ T curFlow; a=b,1 : 0;} vector<pair<int,int>> &edges){ // N+MsqrtN
auto &edge = circ.F.g[i.s+2][cnt[i.s]]; template<typename T>vector<T>Hungarian(const int m = edges.size();
if(i.s != i.e) curFlow = i.r - (edge.c - edge.f); vector<vector<T>>&C){ vector<int> deg(n), pos(n), ord; ord.reserve(n);
else curFlow = i.r; const int J = C.size(), W = C[0].size(); assert(J <= W); vector<vector<int>> gph(n), que(m+1), vec(n);
int srev = sz(g[i.e]), erev = sz(g[i.s]); vector<int> job(W+1, -1); //job[i] - i(worker) matched vector<vector<tuple<int,int,int>>> tri(n);
if(i.s == i.e) srev++; vector<T> ys(J), yt(W + 1), answers;//W-th worker is dummy vector<tuple<int,int,int>> res;
g[i.s].push_back(i.e,srev,+i.r-curFlow,+i.r,+i.cost*(n+1)); const T inf = numeric_limits<T>::max(); for(auto [u,v] : edges) deg[u]++, deg[v]++;
g[i.e].push_back(i.s,erev,-i.l+curFlow,-i.l,-i.cost*(n+1)); for(int j_cur=0; j_cur<J; j_cur++){ for(int i=0; i<n; i++) que[deg[i]].push_back(i);
eps = max(eps, abs(i.cost) * (n+1)); int w_cur = W; job[w_cur] = j_cur; for(int i=m; i>=0; i--) ord.insert(ord.end(), que[i].begin(),
if(i.s != i.e) cnt[i.s] += 2, cnt[i.e] += 2; vector<T> min_to(W+1,inf);vector<int> prv(W+1, -1),in(W+1); que[i].end());
} while(job[w_cur] != -1){ for(int i=0; i<n; i++) pos[ord[i]] = i;
while(true){ eps=0; auto cost=[&](Edge &e, int s, int t){ in[w_cur]=1; T delta=inf; int j = job[w_cur], w_next; for(auto [u,v] : edges) gph[pos[u]].push_back(pos[v]),
return e.cost + p[s] - p[t]; }; for(int w=0; w<W; w++){ if(in[w] != 0) continue; gph[pos[v]].push_back(pos[u]);
for(int i = 0; i < n; i++) for(auto &e : g[i]) if(ckmin(min_to[w], C[j][w]-ys[j]-yt[w])) prv[w]=w_cur; for(int i=0; i<n; i++){
if(e.rem > 0) eps = max(eps, -cost(e, i, e.pos)); if(ckmin(delta, min_to[w])) w_next = w; for(auto j : gph[i]){
if(eps <= T(1)) break; } if(i > j) continue;
eps = max(T(1), eps >> SCALE); for(int w=0; w<=W; w++){ for(int x=0, y=0; x<vec[i].size() && y<vec[j].size(); ){
vector<T> excess(n); queue<int> que; if(in[w] == 0) min_to[w] -= delta; if(vec[i][x] == vec[j][y]) res.emplace_back(ord[i],
auto push = [&](Edge &e, int src, T flow){ else ys[job[w]] += delta, yt[w] -= delta; ord[j], ord[vec[i][x]]), x++, y++;
e.rem -= flow; g[e.pos][e.rev].rem += flow; } /*end for w*/ w_cur = w_next; } /* end while */ else if(vec[i][x] < vec[j][y]) x++; else y++;
excess[src] -= flow; excess[e.pos] += flow; for(int w; w_cur!=-1; w_cur=w)job[w_cur]=job[w=prv[w_cur]]; }
if(excess[e.pos] <= flow && excess[e.pos] > 0) answers.push_back(-yt[W]); vec[j].push_back(i);
que.push(e.pos); } return answers; } }
}; vector<int> ptr(n); }
auto relabel = [&](int v){ 3.17 O(V 3 ) Global Min Cut for(auto &[u,v,w] : res){
ptr[v] = 0; p[v] = -INF; template<typename T, T INF>// 0-based, adj matrix if(pos[u] < pos[v]) swap(u, v);
for(auto &e : g[v]) pair<T, vector<int>> GetMinCut(vector<vector<T>> g){ if(pos[u] < pos[w]) swap(u, w);
if(e.rem>0) p[v] = max(p[v], p[e.pos]-e.cost-eps); int n=g.size(); vector<int> use(n), cut, mn_cut; T mn=INF; if(pos[v] < pos[w]) swap(v, w);
}; for(int phase=n-1; phase>=0; phase--){ tri[u].emplace_back(u, v, w);
for(int i = 0; i < n; i++) for(auto &j : g[i]) vector<int> w=g[0], add=use; int k=0, prv; }
if(j.rem>0 && cost(j, i, j.pos)<0) push(j, i, j.rem); for(int i=0; i<phase; i++){ prv = k; k = -1; res.clear();
while(sz(que)){ for(int j=1; j<n; j++) if(!add[j] && (k==-1 || w[j] > for(int i=n-1; i>=0; i--) res.insert(res.end(),
int x = que.front(); que.pop(); w[k])) k=j; tri[ord[i]].begin(), tri[ord[i]].end());
while(excess[x] > 0){ if(i + 1 < phase){ return res;
for(; ptr[x] < sz(g[x]); ptr[x]++){ for(int j=0; j<n; j++) w[j] += g[k][j]; }
Edge &e = g[x][ptr[x]]; add[k] = 1; continue; } bitset<500> B[500]; // N3/w
if(e.rem > 0 && cost(e, x, e.pos) < 0){ for(int j=0; j<n; j++) g[prv][j] += g[k][j]; long long Count3Cycle(int n, const vector<pair<int,int>>
push(e, x, min(e.rem, excess[x])); for(int j=0; j<n; j++) g[j][prv] = g[prv][j]; &edges){
if(excess[x] == 0) break; use[k] = 1; cut.push_back(k); long long res = 0;
} /* if end */ } /* for end*/ if(w[k] < mn) mn_cut = cut, mn = w[k]; for(int i=0; i<n; i++) B[i].reset();
Soongsil University – PS akgwi Page 12 of 25
for(auto [u,v] : edges) B[u].set(v), B[v].set(u); res.push_back({x, i, j}); vector<int> r(n), val(n), idom(n, -1), sdom(n, -1), o, p(n),
for(int i=0; i<n; i++) for(int j=i+1; j<n; j++) return res; } u(n);
if(B[i].test(j)) res += (B[i] & B[j]).count(); iota(all(r), 0); iota(all(val), 0);
return res / 3; 3.21 O(V E) Shortest Mean Cycle for(int i=0; i<n; i++) for(auto j : g[i]) rg[j].push_back(i);
} template<typename T, T INF> vector<int> // T = V*E*max(C) function<int(int)> find = [&](int v){
// O(n + m * sqrt(m) + th) for graphs without loops or min_mean_cycle(int n, const vector<tuple<int,int,T>> &edges){ if(v == r[v]) return v;
multiedges vector<vector<T>>dp(n+1,vector<T>(n,INF)); // int support! int ret = find(r[v]);
void Find4Cycle(int n, const vector<array<int, 2>> &edge, auto vector<vector<int>>pe(n+1,vector<int>(n,-1)); if(sdom[val[v]] > sdom[val[r[v]]]) val[v] = val[r[v]];
process, int th = 1){ fill(dp[0].begin(),dp[0].end(),0); //0-based,directed return r[v] = ret;
int m = (int)edge.size(); for(int x=1; x<=n; x++){ int id=0; // bellman };
vector<int> deg(n), order, pos(n); for(auto [u,v,w] : edges){ function<void(int)> dfs = [&](int v){
vector<vector<int>> appear(m+1), adj(n), found(n); if(dp[x-1][u] != INF && dp[x-1][u] + w < dp[x][v]) sdom[v] = o.size(); o.push_back(v);
for(auto [u, v]: edge) ++deg[u], ++deg[v]; dp[x][v] = dp[x-1][u] + w, pe[x][v] = id; for(auto i : g[v]) if(sdom[i] == -1) p[i] = v, dfs(i);
for(auto u=0; u<n; u++) appear[deg[u]].push_back(u); id++; } // range based for end! };
for(auto d=m; d>=0; d--) order.insert(order.end(), } T p=1; int q=0, src=-1; //fraction dfs(src); reverse(all(o));
appear[d].begin(), appear[d].end()); for(auto u=0; u<n; u++){ if(dp[n][u] == INF) continue; for(auto &i : o){
for(auto i=0; i<n; i++) pos[order[i]] = i; T cp=-1, cq=0; // ↓ overflow!!! if(sdom[i] == -1) continue;
for(auto i=0; i<m; i++){ for(int x=0;x<=n;x++)if(cp*(n-x) < (dp[n][u]-dp[x][u])*cq) for(auto j : rg[i]){
int u = pos[edge[i][0]], v = pos[edge[i][1]]; cp = dp[n][u] - dp[x][u], cq = n - x; if(sdom[j] == -1) continue;
adj[u].push_back(v), adj[v].push_back(u); if(p * cq > cp * q) src = u, p = cp, q = cq; int x = val[find(j), j];
} } if(src == -1) return {}; if(sdom[i] > sdom[x]) sdom[i] = sdom[x];
T res = 0; vector<int> cnt(n); vector<int> res, po(n, -1); }
for(auto u=0; u<n; u++){ for(int u=src, x = n; ; u=get<0>(edges[pe[x--][u]])){ buf[o[o.size() - sdom[i] - 1]].push_back(i);
for(auto v: adj[u]) if(u < v) for(auto w: adj[v]) if(u < w) if(po[u] != -1)return for(auto j : buf[p[i]]) u[j] = val[find(j), j];
cnt[w] = 0; vector<int>{res.rbegin(),res.rend()-po[u]}; buf[p[i]].clear();
for(auto v: adj[u]) if(u < v) for(auto w: adj[v]) if(u < w) po[u] = res.size(); res.push_back(pe[x][u]); r[i] = p[i];
res += cnt[w] ++; } assert(false); }
} } // return edge index reverse(all(o)); idom[src] = src;
for(auto u=0; u<n; u++){ for(auto i : o){ // WARNING : if different, takes idom
for(auto v: adj[u]) if(u < v) for(auto w: adj[v]) if(u < w) 3.22 O(V 2 ) Stable Marriage Problem if(i != src) idom[i] = sdom[i] == sdom[u[i]] ? sdom[i] :
found[w].clear(); idom[u[i]];
// man : 1~n, woman : n+1~2n
for(auto v: adj[u]) if(u < v) for(auto w: adj[v]) if(u < w) }
struct StableMarriage{
{ for(auto i : o) if(i != src) idom[i] = o[idom[i]];
int n; vector<vector<int>> g;
for(auto x: found[w]){ return idom; // unreachable -> ret[i] = -1
StableMarriage(int n) : n(n), g(2*n+1) { for(int i=1; i<=n+n;
if(!th--) return; }
i++) g[i].reserve(n); }
process(order[u], order[v], order[w], order[x]);
void addEdge(int u, int v){ g[u].push_back(v); } // insert in
}
decreasing order of preference.
found[w].push_back(v);
vector<int> run(){
}
queue<int> q; vector<int> match(2*n+1), ptr(2*n+1);
3.24 O(V E) Vizing Theorem
}
for(int i=1; i<=n; i++) q.push(i); // Graph coloring with (max-degree)+1 colors, O(N^2)
}
while(q.size()){ int C[MX][MX] = {}, G[MX][MX] = {}; // MX ~= 2500
3.20 O(V log V ) Rectlinear MST int i = q.front(); q.pop(); void solve(vector<pii> &E, int N, int M){
template<class T> vector<tuple<T, int, int>> for(int &p=ptr[i]; p<g[i].size(); p++){ int X[MX] = {}, a, b;
rectilinear_minimum_spanning_tree(vector<point<T>> a){ int j = g[i][p]; auto update = [&](int u){ for(X[u] = 1; C[u][X[u]]; X[u]++);
int n = a.size(); vector<int> ind(n); if(!match[j]){ match[i] = j; match[j] = i; break; } };
iota(ind.begin(),ind.end(),0); vector<tuple<T,int,int>> edge; int m = match[j], u = -1, v = -1; auto color = [&](int u, int v, int c){
for(int k=0; k<4; k++){ map<T, int> mp; for(int k=0; k<g[j].size(); k++){ int p = G[u][v]; G[u][v] = G[v][u] = c;
sort(ind.begin(), ind.end(), [&](int i,int j){ if(g[j][k] == i) u = k; if(g[j][k] == m) v = k; C[u][c] = v; C[v][c] = u; C[u][p] = C[v][p] = 0;
return a[i].x-a[j].x < a[j].y-a[i].y;}); } if( p ) X[u] = X[v] = p; else update(u), update(v);
for(auto i: ind){ if(u < v){ return p; }; // end of function : color
for(auto it=mp.lower_bound(-a[i].y); it!=mp.end(); match[m] = 0; q.push(m); match[i] = j; match[j] = i; auto flip = [&](int u, int c1, int c2){
it=mp.erase(it)){ break; int p = C[u][c1], q = C[u][c2];
int j = it->second; point<T> d = a[i] - a[j]; } /*if u < v*/ } /*for-p*/ } /*while*/ swap(C[u][c1], C[u][c2]);
if(d.y > d.x) break; edge.push_back({d.x+d.y,i,j}); return match; } /*vector<int> run*/ if( p ) G[u][p] = G[p][u] = c2;
} }; if( !C[u][c1] ) X[u] = c1; if( !C[u][c2] ) X[u] = c2;
mp.insert({-a[i].y, i}); return p; }; // end of function : flip
} 3.23 O((V + E) log V ) Dominator Tree for(int i = 1; i <= N; i++) X[i] = 1;
for(auto &p: a) if(k & 1) p.x = -p.x; else swap(p.x, p.y); vector<int> DominatorTree(const vector<vector<int>> &g, int for(int t = 0; t < E.size(); t++){
} /*for-k end*/ sort(edge.begin(), edge.end()); src){ // // 0-based int u=E[t].first, v0=E[t].second, v=v0, c0=X[u], c=c0, d;
disjoint_set dsu(n); vector<tuple<T, int, int>> res; int n = g.size(); vector<pii> L; int vst[MX] = {};
for(auto [x, i, j]: edge) if(dsu.merge(i, j)) vector<vector<int>> rg(n), buf(n); while(!G[u][v0]){
Soongsil University – PS akgwi Page 13 of 25
L.emplace_back(v, d = X[v]); for(int i=0; i<e.size(); i++) push(h[e[i].v], {e[i].w, i}); heap = vector<node*>(n, nullptr);
if(!C[v][c]) for(a = (int)L.size()-1; a >= 0; a--) c = vector<int> a(n*2), v(n*2, -1), pa(n*2, -1), r(n*2); priority_queue<pair<ll,ll>,vector<pair<ll,ll>>,greater<>>pq;
color(u, L[a].first, c); iota(a.begin(), a.end(), 0); auto enqueue = [&](int v, ll c, int pa, int pe){
else if(!C[u][d])for(a=(int)L.size()-1;a>=0;a--) auto o = [&](int x) { int y; for(y=x; a[y]!=y; y=a[y]);;; if(dist[v] > c) dist[v] = c, par[v] = pa, pae[v] = pe,
color(u,L[a].first,L[a].second); for(int ox=x; x!=y; ox=x) x = a[x], a[ox] = y; pq.emplace(c, v);
else if( vst[d] ) break; return y; }; }; enqueue(snk, 0, -1, -1); vector<int> ord;
else vst[d] = 1, v = C[u][d]; v[root] = n + 1; int pc = n; while(!pq.empty()){
} for(int i=0; i<n; i++) if(v[i] == -1) { auto [c,v] = pq.top(); pq.pop(); if(dist[v] != c) continue;
if( !G[u][v0] ){ for(int p=i; v[p]==-1 || v[p]==i; p=o(e[r[p]].u)){ ord.push_back(v); for(auto e : rev[v]) enqueue(e.v, c+e.c,
for(;v; v = flip(v, c, d), swap(c, d)); if(v[p] == i){ int q = p; p = pc++; v, e.i);
if(C[u][c0]){ do{ h[q].second = -h[q].first.top().first; }
for(a=(int)L.size()-2; a>=0 && L[a].second!=c; a--); join(h[pa[q]=a[q]=p], h[q]); for(auto &v : ord){
for(; a >= 0; a--) color(u, L[a].first, L[a].second); }while((q=o(e[r[q]].u)) != p); if(par[v] != -1) heap[v] = copy(heap[par[v]]);
} else t--; } v[p] = i; for(auto &e : gph[v]){
} while(!h[p].first.empty() && o(e[top(h[p]).second].u) == if(e.i == pae[v]) continue;
} p) h[p].first.pop(); ll delay = dist[e.v] + e.c - dist[v];
} r[p] = top(h[p]).second; if(delay < 1e18) heap[v] = merge(heap[v], new
} node(make_pair(delay, e.v)));
3.25 O(E + V 3 + V 3T + V 2 2T ) Minimum Steiner Tree } }
vector<int> ans; }
struct SteinerTree{ // O(E + V^3 + V 3^T + V^2 2^T)
for(int i=pc-1; i>=0; i--) if(i != root && v[i] != n) { }
constexpr static int V = 33, T = 8;
for(int f=e[r[i]].v; f!=-1 && v[f]!=n; f=pa[f]) v[f] = n;;; vector<ll> run(int s, int e, int k){
int n, G[V][V], D[1<<T][V], tmp[V];
ans.push_back(r[i]); using state = pair<ll, node*>; dijkstra(e); vector<ll> ans;
void init(int _n){ n = _n;
} priority_queue<state, vector<state>, greater<state>> pq;
memset(G, 0x3f, sizeof G); for(int i=0; i<n; i++) G[i][i]=0;
return ans; if(dist[s] > 1e18) return vector<ll>(k, -1);
} void shortest_path(){ /*floyd 0..n-1*/ }
} ans.push_back(dist[s]);
void add_edge(int u, int v, int w){
if(heap[s]) pq.emplace(dist[s] + heap[s]->val.first, heap[s]);
G[u][v] = G[v][u] = min(G[v][u], w); } 3.27 O(E log V + K log K) K Shortest Walk while(!pq.empty() && ans.size() < k){
int solve(const vector<int>& ter){ int rnd(int l, int r){ /* return random int [l,r] */ } auto [cst, ptr] = pq.top(); pq.pop(); ans.push_back(cst);
int t = (int)ter.size(); memset(D, 0x3f, sizeof D); struct node{ // weight>=0, allow multi edge, self loop for(int j=0; j<2; j++) if(ptr->son[j])
for(int i=0; i<n; i++) D[0][i] = 0; array<node*, 2> son; pair<ll, ll> val; pq.emplace(cst-ptr->val.first + ptr->son[j]->val.first,
for(int msk=1; msk<(1 << t); msk++){ node() : node(make_pair(-1e18, -1e18)) {} ptr->son[j]);
if(msk == (msk & (-msk))){ int who = __lg(msk); node(pair<ll, ll> val) : node(nullptr, nullptr, val) {} int v = ptr->val.second;
for(int i=0; i<n; i++) D[msk][i] = G[ter[who]][i]; node(node *l, node *r, pair<ll,ll> val):son({l,r}),val(val){} if(heap[v]) pq.emplace(cst + heap[v]->val.first, heap[v]);
continue; }; }
} node* copy(node *x){ return x ? new node(x->son[0], x->son[1], while(ans.size() < k) ans.push_back(-1);
for(int i=0; i<n; i++) x->val) : nullptr; } return ans;
for(int sub=(msk-1)&msk; sub; sub=(sub-1)&msk) node* merge(node *x, node *y){ // precondition: x, y both points }
D[msk][i] = min(D[msk][i], D[sub][i] + D[msk^sub][i]); to new entity
memset(tmp, 0x3f, sizeof tmp); if(!x || !y) return x ? x : y;
for(int i=0; i<n; i++) for(int j=0; j<n; j++) if(x->val > y->val) swap(x, y);
tmp[i] = min(tmp[i], D[msk][j] + G[j][i]); int rd = rnd(0, 1); if(x->son[rd])
3.28 O(V + E) Chordal Graph, Tree Decomposition
for(int i=0; i<n; i++) D[msk][i] = tmp[i]; x->son[rd]=copy(x->son[rd]); struct Set { list<int> L; int last; Set() { last = 0; } };
} x->son[rd] = merge(x->son[rd], y); return x; struct PEO {
return *min_element(D[(1<<t)-1], D[(1<<t)-1]+n); } int N; list<Set> L;
} struct edge{ vector<vector<int>> g; vector<int> vis, res;
}; ll v, c, i; edge() = default; vector<list<Set>::iterator> ptr;
edge(ll v, ll c, ll i) : v(v), c(c), i(i) {} vector<list<int>::iterator> ptr2;
3.26 O(E log V ) Directed MST }; PEO(int n, vector<vector<int> > _g) {
using D = int; struct edge { int u, v; D w; }; vector<vector<edge>> gph, rev; int idx; N = n; g = _g;
vector<int> DirectedMST(vector<edge> &e, int n, int root){ void init(int n){ gph = rev = vector<vector<edge>>(n); idx=0; } for (int i = 1; i <= N; i++) sort(g[i].begin(), g[i].end());
using T = pair<D, int>; // 0-based, return index of edges void add_edge(int s, int e, ll x){ vis.resize(N + 1); ptr.resize(N + 1); ptr2.resize(N + 1);
using PQ = pair<priority_queue <T,vector<T>,greater<T>>, D>; gph[s].emplace_back(e, x, idx); L.push_back(Set());
auto push = [](PQ &pq, T v){ rev[e].emplace_back(s, x, idx); for (int i = 1; i <= N; i++) {
pq.first.emplace(v.first-pq.second, v.second); }; assert(x >= 0); idx++; L.back().L.push_back(i);
auto top = [](const PQ &pq) -> T { } ptr[i] = L.begin(); ptr2[i] = prev(L.back().L.end());
auto r = pq.first.top(); return {r.first + pq.second, vector<int> par, pae; vector<ll> dist; vector<node*> heap; }
r.second}; }; void dijkstra(int snk){ // replace this to SPFA if edge weight }
auto join = [&push, &top](PQ &a, PQ &b) { is negative pair<bool, vector<int>> Run() {
if(a.first.size() < b.first.size()) swap(a, b); int n = gph.size(); // lexicographic BFS
while(!b.first.empty()) push(a, top(b)), b.first.pop(); }; par = pae = vector<int>(n, -1); int time = 0;
vector<PQ> h(n * 2); dist = vector<ll>(n, 0x3f3f3f3f3f3f3f3f); while (!L.empty()) {
Soongsil University – PS akgwi Page 14 of 25
if (L.front().L.empty()) { L.pop_front(); continue; } for(int next : g[n]) if (vis[next] && mn > ord[next]) mn = void upd(int u, int v){ if (!sl[v] || d(e[u][v]) <
auto it = L.begin(); ord[next], idx = next; d(e[sl[v]][v])) sl[v] = u; }
int n = it->L.front(); it->L.pop_front(); assert(idx != -1); idx = p[idx]; void ss(int v){
vis[n] = ++time; // 두 set인 V[idx]와 g[n](visited ver)가 같나? sl[v]=0; for(int u=1; u<=n; u++) if(e[u][v].w > 0 && st[u]
res.push_back(n); // V[idx]의 모든 원소가 g[n]에서 나타나는지 판별로 충분하다. != v && !s[st[u]]) upd(u, v);
for (int next : g[n]) { int die = 0; }
if (vis[next]) continue; for(int x : V[idx]) { void ins(int u){ if(u <= n) q[++t] = u; else for(int v : p[u])
if (ptr[next]->last != time) { if (!binary_search(g[n].begin(), g[n].end(), x)) { die = ins(v); }
L.insert(ptr[next], Set()); ptr[next]->last = time; 1; break; } void mdf(int u, int w){ st[u]=w; if(u > n) for(int v : p[u])
} } mdf(v, w); }
ptr[next]->L.erase(ptr2[next]); ptr[next]--; if (!die) { V[idx].push_back(n), p[n] = idx; } // 기존 집합에 int gr(int u,int v){
ptr[next]->L.push_back(next); 추가 if ((v=find(p[u].begin(), p[u].end(), v) - p[u].begin()) &
ptr2[next] = prev(ptr[next]->L.end()); else { // 새로운 집합을 자식으로 추가 1){
} P++; reverse(p[u].begin()+1, p[u].end()); return
} G[idx].push_back(P); // 자식으로만 단방향으로 잇자. (int)p[u].size() - v;
// PEO existence check V[P].push_back(n); }
for (int n = 1; n <= N; n++) { for(int next : g[n]) if (vis[next]) V[P].push_back(next); return v; }
int mx = 0; p[n] = P; void stm(int u, int v){
for (int next : g[n]) if (vis[n] > vis[next]) mx = max(mx, } lk[u] = e[u][v].v;
vis[next]); } if(u <= n) return; Q w = e[u][v];
if (mx == 0) continue; for(int i=1; i<=P; i++) sort(V[i].begin(), V[i].end()); int x = b[u][w.u], y = gr(u,x);
int w = res[mx - 1]; } for(int i=0; i<y; i++) stm(p[u][i], p[u][i^1]);
for (int next : g[n]) { stm(x,v);rotate(p[u].begin(), p[u].begin()+y, p[u].end()); }
if (vis[w] > vis[next] && !binary_search(g[w].begin(), 3.29 O(V 3 ) General Matching void aug(int u, int v){
g[w].end(), next)){ int N, M, R, Match[555], Par[555], Chk[555], Prv[555], Vis[555]; int w = st[lk[u]]; stm(u, v); if (!w) return;
vector<int> chk(N+1), par(N+1, -1); // w와 next가 vector<int> G[555]; // n 500 20ms stm(w, st[f[w]]); aug(st[f[w]], w); }
이어져 있지 않다면 not chordal int Find(int x){return x == Par[x] ? x : Par[x] = Find(Par[x]);} int lca(int u, int v){
deque<int> dq{next}; chk[next] = 1; int LCA(int u, int v){ static int cnt = 0; for(++id; u|v; swap(u, v)){
while (!dq.empty()) { for(cnt++; Vis[u]!=cnt; swap(u, v)) if(u) Vis[u] = cnt, u = if(!u) continue; if(ed[u] == id) return u;
int x = dq.front(); dq.pop_front(); Find(Prv[Match[u]]); ed[u] = id; if(u = st[lk[u]]) u = st[f[u]]; // not ==
for (auto y : g[x]) { return u; } }
if (chk[y] || y == n || y != w && void Blossom(int u, int v, int rt, queue<int> &q){ return 0; }
binary_search(g[n].begin(), g[n].end(), y)) for(; Find(u)!=rt; u=Prv[v]){ void add(int u, int a, int v){
continue; Prv[u] = v; Par[u] = Par[v=Match[u]] = rt; int x = n+1; while(x <= m && st[x]) x++;
dq.push_back(y); chk[y] = 1; par[y] = x; if(Chk[v] & 1) q.push(v), Chk[v] = 2; if(x > m) m++;
} } } lab[x] = s[x] = st[x] = 0; lk[x] = lk[a];
} bool Augment(int u){ // iota Par 0, fill Chk 0 p[x].clear(); p[x].push_back(a);
vector<int> cycle{next, n}; queue<int> Q; Q.push(u); Chk[u] = 2; for(int i=u, j; i!=a; i=st[f[j]]) p[x].push_back(i),
for (int x=w; x!=next; x=par[x]) cycle.push_back(x); while(!Q.empty()){ u = Q.front(); Q.pop(); p[x].push_back(j=st[lk[i]]), ins(j);
return {false, cycle}; for(auto v : G[u]){ reverse(p[x].begin()+1, p[x].end());
} if(Chk[v] == 0){ for(int i=v, j; i!=a; i=st[f[j]]) p[x].push_back(i),
} Prv[v]=u; Chk[v]=1; Q.push(Match[v]); Chk[Match[v]]=2; p[x].push_back(j=st[lk[i]]), ins(j);
} if(!Match[v]){ for(; u; v=u) u = Match[Prv[v]], mdf(x,x); for(int i=1; i<=m; i++) e[x][i].w=e[i][x].w=0;
reverse(res.begin(), res.end()); Match[Match[v]=Prv[v]] = v;;; return true; } memset(b[x]+1, 0, n*sizeof b[0][0]);
return {true, res}; } for (int u : p[x]){
} else if(Chk[v] == 2){ int l = LCA(u, v); Blossom(u, v, l, for(v=1; v<=m; v++) if(!e[x][v].w || d(e[u][v]) <
}; Q), Blossom(v, u, l, Q); } d(e[x][v])) e[x][v] = e[u][v],e[v][x] = e[v][u];
bool vis[200201]; // 배열 크기 알아서 수정하자. } /* for v */ } /* while */ for(v=1; v<=n; v++) if(b[u][v]) b[x][v] = u;
int p[200201], ord[200201], P = 0; // P=정점 개수 return 0; } }
vector<int> V[200201], G[200201]; // V=bags, G=edges void Run(){ for(int i=1; i<=N; i++) if(!Match[i]) R += ss(x); }
void tree_decomposition(int N, vector<vector<int> > g) { Augment(i); } void ex(int u){ // s[u] == 1
for(int i=1; i<=N; i++) sort(g[i].begin(), g[i].end()); for(int x : p[u]) mdf(x, x);
vector<int> peo = PEO(N, g).Run(), rpeo = peo; 3.30 O(V 3 ) Weighted General Matching int a = b[u][e[u][f[u]].u],r = gr(u, a);
reverse(rpeo.begin(), rpeo.end()); namespace weighted_blossom_tree{ // n 400 w 1e8 700ms, n 500 w for(int i=0; i<r; i+=2){
for(int i=0; i<peo.size(); i++) ord[peo[i]] = i; 1e6 300ms int x = p[u][i], y = p[u][i+1];
for(int n : rpeo) { // tree decomposition #define d(x) (lab[x.u]+lab[x.v]-e[x.u][x.v].w*2) f[x] = e[y][x].u; s[x] = 1; s[y] = 0; sl[x] = 0; ss(y);;;
vis[n] = true; const int N=403*2; using ll = long long; using T = int; // sum ins(y); }
if (n == rpeo[0]) { // 처음 of weight, single weight s[a] = 1; f[a] = f[u];
P++; V[P].push_back(n); p[n] = P; continue; const T inf=numeric_limits<T>::max()>>1; for(int i=r+1;i<p[u].size();i++)s[p[u][i]]=-1, ss(p[u][i]);
} struct Q{ int u, v; T w; } e[N][N]; vector<int> p[N]; st[u] = 0; }
int mn = INF, idx = -1; int n, m=0, id, h, t, lk[N], sl[N], st[N], f[N], b[N][N],
s[N], ed[N], q[N]; T lab[N];
Soongsil University – PS akgwi Page 15 of 25
bool on(const Q &e){ 4 Math };
int u=st[e.u], v=st[e.v], a;
4.1 Binary GCD, Extend GCD, CRT, Combination
if(s[v] == -1) f[v] = e.u, s[v] = 1, a = st[lk[v]], sl[v] = 4.2 Partition Number
sl[a] = s[a] = 0, ins(a); ll binary_gcd(ll a, ll b){
for(int j=1; j*(3*j-1)/2<=i; j++) P[i] +=
else if(!s[v]){ if(a == 0 || b == 0) return a + b;
(j%2?1:-1)*P[i-j*(3*j-1)/2], P[i] %= MOD;
a = lca(u, v); if(!a) return aug(u,v), aug(v,u), true; int az = __builtin_ctzll(a), bz = __builtin_ctzll(b);
for(int j=1; j*(3*j+1)/2<=i; j++) P[i] +=
else add(u,a,v); int shift = min(az, bz); b >>= bz;
(j%2?1:-1)*P[i-j*(3*j+1)/2], P[i] %= MOD;
} while(a){ a >>= az; ll diff = b-a;
vector<ModInt> res(sz+1); res[0] = 1; int sq=sqrt(sz);
return false; } az = __builtin_ctz(diff); b = min(a, b); a = abs(diff);
vector<vector<ModInt>> p(2, vector<ModInt>(sz+1)), d=p;
bool bfs(){ } return b << shift;
for(int k=1; k<sq; k++){ p[0][0] = k == 1; // calc p[k][n]
memset(s+1, -1, m*sizeof s[0]); memset(sl+1, 0, m*sizeof } // return [g,x,y] s.t. ax+by=gcd(a,b)=g
for(int n=1; n<=sz; n++){
sl[0]); tuple<ll,ll,ll> ext_gcd(ll a, ll b){
p[k&1][n] = p[k-1&1][n-1] + (n-k>=0 ? p[k&1][n-k] : 0);
h = 1; t = 0; for(int i=1; i<=m; i++) if(st[i] == i && if(b == 0) return {a, 1, 0}; auto [g,x,y] = ext_gcd(b, a % b);
res[n] += p[k&1][n]; }}
!lk[i]) f[i] = s[i] = 0, ins(i); return {g, y, x - a/b * y}; }
for(int a=sq; a>=0; a--) for(int b=sq; b<=sz; b++){
if(h > t) return 0; ll inv(ll a, ll m){ //return x when ax mod m = 1, fail -> -1
d[a&1][b] = d[a+1&1][b-sq] + p[sq-1&1][b-1] + (b-a-1>=0 ?
while (true){ auto [g,x,y] = ext_gcd(a, m); return g == 1 ? mod(x, m) : -1;}
d[a&1][b-a-1] : 0);
while (h <= t){ void DivList(ll n){ // {n/1, n/2, ... , n/n}, size <= 2 sqrt n
if(a == 0) res[b] += d[a&1][b];
int u = q[h++]; for(ll i=1, j=1; i<=n; i=j+1) Report(i, j=n/(n/i), n/i); }
}
if (s[st[u]] != 1) for (int v=1; v<=n; v++) if void Div2List(ll n){// n/(i^2), n^{3/4}
(e[u][v].w > 0 && st[u] != st[v]) for(ll i=1, j=1; i*i<=n; i=j+1){
j = (ll)floorl(sqrtl(n/(n/(i*i)))); Report(i, j, n/(i*i));
4.3 Diophantine
if(d(e[u][v])) upd(u, st[v]); else if(on(e[u][v]))
} }//square free: sum_{i=1..sqrt n} mu(i)floor(n/(i^2)) // solutions to ax + by = c where x in [xlow, xhigh] and y in
return true;
pair<ll,ll> crt(ll a1, ll m1, ll a2, ll m2){ [ylow, yhigh]
}
ll g = gcd(m1, m2), m = m1 / g * m2; // cnt, leftsol, rightsol, gcd of a and b
T x = inf;
if((a2 - a1) % g) return {-1, -1}; template<class T> array<T, 6> solve_linear_diophantine(T a, T b,
for(int i=n+1; i<=m; i++) if(st[i] == i && s[i] == 1) x =
ll md = m2/g, s = mod((a2-a1)/g, m2/g); T c, T xlow, T xhigh, T ylow, T yhigh){
min(x, lab[i]>>1);
ll t = mod(get<1>(ext_gcd(m1/g%md, m2/g)), md); T g, x, y = euclid(a >= 0 ? a : -a, b >= 0 ? b : -b, x, y);
for(int i=1; i<=m; i++) if(st[i] == i && sl[i] && s[i] !=
return { a1 + s * t % md * m1, m }; } array<T, 6> no_sol{0, 0, 0, 0, 0, g};
1) x = min(x, d(e[sl[i]][i])>>s[i]+1);
pair<ll,ll> crt(const vector<ll> &a, const vector<ll> &m){ if(c % g) return no_sol; x *= c / g, y *= c / g;
for(int i=1; i<=n; i++) if(~s[st[i]]) if((lab[i] +=
ll ra = a[0], rm = m[0]; if(a < 0) x = -x; if(b < 0) y = -y;
(s[st[i]]*2-1)*x) <= 0) return false;
for(int i=1; i<m.size(); i++){ a /= g, b /= g, c /= g;
for(int i=n+1 ;i<=m; i++) if(st[i] == i && ~s[st[i]])
auto [aa,mm] = crt(ra, rm, a[i], m[i]); auto shift = [&](T &x, T &y, T a, T b, T cnt){ x += cnt * b,
lab[i] += (2-s[st[i]]*4)*x;
if(mm == -1) return {-1, -1}; else tie(ra,rm) = tie(aa,mm); y -= cnt * a; };
h = 1; t = 0;
} return {ra, rm}; } int sign_a = a > 0 ? 1 : -1, sign_b = b > 0 ? 1 : -1;
for(int i=1; i<=m; i++) if(st[i] == i && sl[i] &&
struct Lucas{ // init : O(P), query : O(log P) shift(x, y, a, b, (xlow - x) / b);
st[sl[i]] != i && !d(e[sl[i]][i]) && on(e[sl[i]][i]))
const size_t P; vector<ll> fac, inv; if(x < xlow) shift(x, y, a, b, sign_b);
return true;
ll Pow(ll a, ll b){ /* return a^b mod P */ } if(x > xhigh) return no_sol;
for(int i=n+1; i<=m; i++) if(st[i] == i && s[i] == 1 &&
Lucas(size_t P):P(P),fac(P),inv(P){ /* init fac, facinv */ } T lx1 = x; shift(x, y, a, b, (xhigh - x) / b);
!lab[i]) ex(i);
ll small(ll n, ll r) const { /* n! / r! / (n-r)! */ } if(x > xhigh) shift(x, y, a, b, -sign_b);
}
ll calc(ll n, ll r) const { if(n<r || n<0 || r<0) return 0; T rx1 = x; shift(x, y, a, b, -(ylow - y) / a);
return 0; }
if(!n || !r || n == r) return 1; if(y < ylow) shift(x, y, a, b, -sign_a);
template<typename TT> pair<int,ll> run(int N, const
else return small(n%P, r%P) * calc(n/P, r/P) % P; } if(y > yhigh) return no_sol;
vector<tuple<int,int,TT>> &edges){ // 1-based
}; T lx2 = x; shift(x, y, a, b, -(yhigh - y) / a);
memset(ed+1, 0, m*sizeof ed[0]); memset(lk+1, 0, m*sizeof
template<ll p, ll e> struct CombinationPrimePower{ if(y > yhigh) shift(x, y, a, b, sign_a);
lk[0]);
vector<ll> val; ll m; // init : O(p^e), query : O(log p) T rx2 = x; if(lx2 > rx2) swap(lx2, rx2);
n = m = N; id = 0; iota(st+1, st+n+1, 1); T wm = 0; ll r =
CombinationPrimePower(){ T lx = max(lx1, lx2), rx = min(rx1, rx2);
0;
m=1; for(int i=0; i<e; i++) m *= p; val.resize(m); val[0]=1; if(lx > rx) return no_sol;
for(int i=1; i<=n; i++) for(int j=1; j<=n; j++) e[i][j] =
for(int i=1; i<m; i++)val[i] = val[i-1] * (i%p ? i : 1) % m; return {(rx - lx) / (b >= 0 ? b : -b) + 1, lx, (c - lx * a)
{i,j,0};
} / b, rx, (c - rx * a) / b, g};
for(auto [u,v,w] : edges) wm = max(wm,
pair<ll,ll> factorial(int n){ if(n < p) return {0, val[n]}; }
e[v][u].w=e[u][v].w=max(e[u][v].w,(T)w));
for(int i=1; i<=n; i++) p[i].clear(); int k = n / p; auto v = factorial(k);
for(int i=1; i<=n; i++) for (int j=1; j<=n; j++) b[i][j] = int cnt = v.first + k, kp = n / m, rp = n % m; 4.4 FloorSum
i*(i==j); ll ret=v.second * Pow(val[m-1], kp%2, m) % m * val[rp] % m; // sum of floor((A*i+B)/M) over 0 <= i < N in O(log(N+M+A+B))
fill_n(lab+1, n, wm); int match = 0; while(bfs()) match++; return {cnt, ret}; } // Also, sum of i * floor((A*i+B)/M) and floor((A*i+B)/M)^2
for(int i=1; i<=n; i++) if(lk[i]) r += e[i][lk[i]].w; ll calc(int n, int r){ if(n < 0 || r < 0 || n < r) return 0; template<class T, class U> // T must be able to hold arg^2
return {match, r/2}; auto v1=factorial(n), v2=factorial(r), v3=factorial(n-r); array<U, 3> weighted_floor_sum(T n, T m, T a, T b){
} ll cnt = v1.first - v2.first - v3.first; array<U, 3> res{}; auto[qa,ra]=div(a,m); auto[qb,rb]=div(b,m);
#undef d ll ret = v1.second * inv(v2.second, m) % m * inv(v3.second, if(T n2 = (ra * n + rb) / m){
} using weighted_blossom_tree::run, weighted_blossom_tree::lk; m) % m; auto prv=weighted_floor_sum<T,U>(n2, ra, m, m-rb-1);
if(cnt >= e) return 0; res[0] += U(n-1)*n2 - prv[0];
for(int i=1; i<=cnt; i++) ret = ret * p % m; res[1] += (U(n-1)*n*n2 - prv[0] - prv[2]) / 2;
return ret; } res[2] += U(n-1)*(n2-1)*n2 - 2*prv[1] + res[0];
Soongsil University – PS akgwi Page 16 of 25
} hi.p += lo.p * adv; hi.q += lo.q * adv; } return {a, rank, det, out}; // linear system: get RREF(A|b)
res[2] += U(n-1)*n*(2*n-1)/6 * qa*qa + U(n)*qb*qb; dir = !dir; swap(lo, hi); A = B; B = adv != 0; } // 0 0 ... 0 b[i]: inconsistent, rank < len(A[0]): multiple
res[2] += U(n-1)*n * qa*qb + 2*res[0]*qb + 2*res[1]*qa; } // get det(A) mod M, M can be composite number
res[0] += U(n-1)*n/2 * qa + U(n)*qb; return dir ? hi : lo; // remove mod M -> get pure det(A) in integer
res[1] += U(n-1)*n*(2*n-1)/6 * qa + U(n-1)*n/2 * qb; } ll Det(vector<vector<ll>> a){//destroy matrix, n500 -400ms
return res; int n = a.size(); ll ans = 1;
} 4.7 O(N 3 log 1/) Polynomial Equation for(int i=0; i<n; i++){
ll modsum(ull to, ll c, ll k, ll m){ for(int j=i+1; j<n; j++){
vector<double> poly_root(vector<double> p, double xmin, double
c = (c % m + m) % m; k = (k % m + m) % m; while(a[j][i] != 0){ // gcd step
xmax){
return to*c + k*sumsq(to) - m*divsum(to, c, k, m); ll t = a[i][i] / a[j][i];
if(p.size() == 2){ return {-p[0] / p[1]}; }
} // sum (ki+c)%m 0<=i<to, O(log m) large constant if(t)for(int k=i;k<n;k++) a[i][k]=(a[i][k]-a[j][k]*t)%M;
vector<double> ret, der(p.size()-1);
swap(a[i], a[j]); ans *= -1;
for(int i=0; i<der.size(); i++) der[i] = p[i+1] * (i + 1);
4.5 XOR Basis (XOR Maximization) auto dr = poly_root(der, xmin, xmax);
}
vector<ll> basis; // ascending }
dr.push_back(xmin-1); dr.push_back(xmax+1);
for(int i=0; i<n; i++){ ll x; cin >> x; ans = ans * a[i][i] % M; if(!ans) return 0;
sort(dr.begin(), dr.end());
for(int j=(int)basis.size()-1; j>=0; j--) x=min(x,basis[j]^x); } return (ans + M) % M;
for(int i=0; i+1<dr.size(); i++){
if(x)basis.insert(lower_bound(basis.begin(),basis.end(), x),x); }
double l = dr[i], h = dr[i+1]; bool sign = calc(p, l) > 0;
}//xor maximization, reverse -> for(auto i:basis)r=max(r,r^i); if (sign ^ (calc(p, h) > 0)){
// minimization, return basis.back(), WARNING: x=0 => return 0 for(int it=0; it<60; it++){ // while(h-l > 1e-8)
// choose 2k element => solve(a1^a2, a1^a3, a1^a4, ...) double m = (l + h) / 2, f = calc(p, m); 4.9 Berlekamp + Kitamasa
if ((f <= 0) ^ sign) l = m; else h = m; const int mod = 1e9+7; ll pw(ll a, ll b){/*a^b mod M*/}
4.6 Stern Brocot Tree } vector<int> berlekamp_massey(vector<int> x){
pair<ll,ll> Solve(ld l, ld r){//find l<p/q<r -> min q -> min p ret.push_back((l + h) / 2); int n = x.size(),L=0,m=0; ll b=1; if(!n) return {};
auto g=[](ll v,pair<ll,ll>a,pair<ll,ll>b)->pair<ll,ll>{ } vector<int> C(n), B(n), T; C[0]=B[0]=1;
return { v * a.first + b.first, v * a.second + b.second }; } for(int i=0; ++m && i<n; i++){ ll d = x[i] % mod;
}; return ret; for(int j=1; j<=L; j++) d = (d + 1LL * C[j] * x[i-j]) % mod;
auto f = [g](ll v, pair<ll,ll> a, pair<ll,ll> b) -> ld { } if(!d) continue; T=C; ll c = d * pw(b, mod-2) % mod;
auto [p,q] = g(v, a, b); return ld(p) / q; }; for(int j=m; j<n; j++) C[j] = (C[j] - c * B[j-m]) % mod;
pair<ll,ll> s(0, 1), e(1, 0); 4.8 Gauss Jordan Elimination if(2 * L <= i) L = i-L+1, B = T, b = d, m = 0;
while(true){ template<typename T> // return {rref, rank, det, inv} }
pair<ll,ll> m(s.first+e.first, s.second+e.second); tuple<vector<vector<T>>, int, T, vector<vector<T>>> C.resize(L+1); C.erase(C.begin());
ld v = 1.L * m.first / m.second; Gauss(vector<vector<T>> a, bool square=true){ // n500 -400ms for(auto &i : C) i = (mod - i) % mod; return C;
if(v >= r){ int n = a.size(), m = a[0].size(), rank = 0;//bitset 4096-700 } // O(NK + N log mod)
ll ks = 1, ke = 1; while(f(ke, s, e) >= r) ke *= 2; vector<vector<T>> out(n, vector<T>(m, 0)); T det = T(1); int get_nth(vector<int> rec, vector<int> dp, ll n){
while(ks <= ke){ for(int i=0; i<n; i++) if(square) out[i][i] = T(1); int m = rec.size(); vector<int> s(m), t(m); ll ret=0;
ll km = (ks + ke) / 2; for(int i=0; i<m; i++){ s[0] = 1; if(m != 1) t[1] = 1; else t[0] = rec[0];
if(f(km, s, e) >= r) ks = km + 1; else ke = km - 1; if(rank == n) break; auto mul = [&rec](vector<int> v, vector<int> w){
} e = g(ke, s, e); if(IsZero(a[rank][i])){ int m = v.size(); vector<int> t(2*m);
} T mx = T(0); int idx = -1; // fucking precision error for(int j=0; j<m; j++) for(int k=0; k<m; k++){
else if(v <= l){ for(int j=rank+1; j<n; j++) if(mx < abs(a[j][i])) mx = t[j+k] = (t[j+k] + 1LL * v[j] * w[k]) % mod;
ll ks = 1, ke = 1; while(f(ke, e, s) <= l) ke *= 2; abs(a[j][i]), idx = j; }
while (ks <= ke){ if(idx == -1 || IsZero(a[idx][i])){ det = 0; continue; } for(int j=2*m-1; j>=m; j--) for(int k=1; k<=m; k++){
ll km = (ks + ke) / 2; for(int k=0; k<m; k++){ t[j-k] = (t[j-k] + 1LL * t[j] * rec[k-1]) % mod;
if(f(km, e, s) <= l) ks = km + 1; else ke = km - 1; a[rank][k] = Add(a[rank][k], a[idx][k]); }
} s = g(ke, e, s); if(square)out[rank][k]=Add(out[rank][k],out[idx][k]); t.resize(m); return t;
} } };
else return m; } for(; n; n>>=1, t=mul(t,t)) if(n & 1) s=mul(s,t);
} det = Mul(det, a[rank][i]); for(int i=0; i<m; i++) ret += 1LL * s[i] * dp[i] % mod;
} T coeff = Div(T(1), a[rank][i]); return ret % mod;
struct Frac { ll p, q; };//find smallest 0 <= p/q <= 1 (p,q<=N) for(int j=0; j<m; j++) a[rank][j] = Mul(a[rank][j], coeff); } // O(N2 log X)
template<class F> Frac fracBS(F f, ll N) { // s.t. f(p/q) true for(int j=0; j<m; j++) if(square) out[rank][j] = int guess_nth_term(vector<int> x, ll n){
bool dir = 1, A = 1, B = 1; // O(log N) Mul(out[rank][j], coeff); if(n < x.size()) return x[n];
Frac lo{0, 1}, hi{1, 1}; // Set hi to 1/0 to search (0, N] for(int j=0; j<n; j++){ vector<int> v = berlekamp_massey(x);
if(f(lo)) return lo; assert(f(hi)); if(rank == j) continue; return v.empty() ? 0 : get_nth(v, x, n);
while(A != 0 || B != 0){ T t = a[j][i]; // Warning: [j][k], [rank][k] }
ll adv = 0, step = 1; // move hi if dir, else lo for(int k=0; k<m; k++) a[j][k] = Sub(a[j][k], struct elem{int x, y, v;}; // A_(x, y) <- v, 0-based. no
for(int si=0; step; (step*=2)>>=si){ adv += step; Mul(a[rank][k], t)); duplicate please..
Frac mid{lo.p * adv + hi.p, lo.q * adv + hi.q}; for(int k=0; k<m; k++) if(square) out[j][k] = vector<int> get_min_poly(int n, vector<elem> M){
if(abs(mid.p)>N || mid.q>N || dir != f(mid)) Sub(out[j][k], Mul(out[rank][k], t)); // smallest poly P such that A^i = sum_{j < i} {A^j \times
adv -= step, si = 2; } P_j}
} rank++; // linear system: warning len(A) != len(A[0]) vector<int> rnd1, rnd2, gobs; mt19937 rng(0x14004);
Soongsil University – PS akgwi Page 17 of 25
auto gen = [&rng](int lb, int ub){ return phi: euler totient function for(ll r=1; ; r++){
uniform_int_distribution<int>(lb, ub)(rng); }; sigma_k: sum of k-th power of divisors bool flag = true; // Warning: 64bit Pow
for(int i=0; i<n; i++) rnd1.push_back(gen(1, mod-1)), sigma = sigma_1, d = tau = sigma_0 for(auto [d,e] : v) if(PowMod(r, (p-1)/d, p) == 1){ flag =
rnd2.push_back(gen(1, mod-1)); sigma_k = id_k * 1 | sigma = id * 1 false; break; }
for(int i=0; i<2*n+2; i++){ int tmp = 0; id_k = sigma_k * mu | id = sigma * mu if(flag) return r;
for(int j=0; j<n; j++) tmp = (tmp + 1LL * rnd2[j] * rnd1[j]) e = 1 * mu | d = 1 * 1 | 1 = d * mu }
% mod; phi * 1 = id | phi = id * mu | sigma = phi * d }
gobs.push_back(tmp); vector<int> nxt(n); g = f * 1 iff f = g * mu */ // Given A, B, P, solve A^x === B mod P, return smallest value
for(auto &j : M) nxt[j.x] = (nxt[j.x] + 1LL * j.v * template<class T, class F1, class F2, class F3> ll DiscreteLog(ll A, ll B, ll P){ // O(sqrt P) with hash set
rnd1[j.y]) % mod; struct xudyh_sieve{ __gnu_pbds::gp_hash_table<ll,__gnu_pbds::null_type> st;
rnd1 = nxt; T th; // threshold, 2(single query) ~ 5 * MAXN^2/3 ll t = ceil(sqrt(P)), k = 1; // use binary search?
} auto v = berlekamp_massey(gobs); F1 pf; F2 pg; F3 pfg; for(int i=0; i<t; i++) st.insert(k), k = k * A % P;
return vector<int>(v.rbegin(), v.rend()); // prefix sum of f(up to th), g(easy to calc), f*g(easy to ll inv = Pow(k, P-2, P);
} calc) for(int i=0, s=1; i<t; i++, s=s*inv%P){
ll det(int n, vector<elem> M){ unordered_map<T, T> mp; // f * g means dirichlet conv. ll x = B * s % P;
vector<int> rnd; mt19937 rng(0x14004); xudyh_sieve(T th,F1 pf,F2 pg,F3 if(st.find(x) == st.end()) continue;
auto gen = [&rng](int lb, int ub){ return pfg):th(th),pf(pf),pg(pg),pfg(pfg){} for(int j=0, f=1; j<t; j++, f=f*A%P){
uniform_int_distribution<int>(lb, ub)(rng); }; // Calculate the preix sum of a multiplicative f up to n if(f == x) return i * t + j;
for(int i=0; i<n; i++) rnd.push_back(gen(1, mod-1)); T query(T n){ // O(n^2/3) }
for(auto &i : M) i.v = 1LL * i.v * rnd[i.y] % mod; if(n <= th) return pf(n); if(mp.count(n)) return mp[n]; }
auto sol = get_min_poly(n, M)[0]; if(n % 2 == 0) sol = mod - T res = pfg(n); return -1;
sol; for(T low = 2, high = 2; low <= n; low = high + 1){ }
for(auto &i : rnd) sol = 1LL * sol * pw(i, mod-2) % mod; high = n / (n / low); // Given A, P, solve X^2 === A mod P, return arbitrary
return sol; res -= (pg(high) - pg(low - 1)) * query(n / low); // MOD ll DiscreteSqrt(ll A, ll P){//O(log^2P), O(logP) in random data
} } if(A == 0) return 0;
return mp[n] = res / pg(1);//Pow(pg(1),MOD-2)? if(Pow(A, (P-1)/2, P) != 1) return -1;
4.10 Linear Sieve } if(P % 4 == 3) return Pow(A, (P+1)/4, P);
// sp : 최소 소인수, 소수라면 0 }; ll s = P - 1, n = 2, r = 0, m;
// tau : 약수 개수, sigma : 약수 합 while(~s & 1) r++, s >>= 1;
// phi : n 이하 자연수 중 n과 서로소인 개수 4.12 Miller Rabin + Pollard Rho while(Pow(n, (P-1)/2, P) != P-1) n++;
// mu : non square free이면 0, 그렇지 않다면 (-1)^(소인수 종류) // 32bit : 2, 7, 61 / ull MulMod, PowMod (cast __uint128_t) ll x = Pow(A, (s+1)/2, P), b = Pow(A, s, P), g = Pow(n, s, P);
// e[i] : 소인수분해에서 i의 지수 // 64bit : 2, 325, 9375, 28178, 450775, 9780504, 1795265022 for(;; r=m){
vector<int> prime; bool MillerRabin(ull n, ull a){ ll t = b; for(m=0; m<r && t!=1; m++) t = t * t % P;
int sp[sz], e[sz], phi[sz], mu[sz], tau[sz], sigma[sz]; if(a % n == 0) return true; int cnt = __builtin_ctzll(n - 1); if(!m) return x;
phi[1] = mu[1] = tau[1] = sigma[1] = 1; ull p=PowMod(a, n>>cnt, n); if(p==1 || p+1==n) return true; ll gs = Pow(g, 1LL << (r-m-1), P);
for(int i=2; i<=n; i++){ while(cnt--) if((p=MulMod(p,p,n)) == n - 1) return true; g = gs * gs % P; x = x * gs % P; b = b * g % P;
if(!sp[i]){ return false; }
prime.push_back(i); } }
e[i] = 1; phi[i] = i-1; mu[i] = -1; tau[i] = 2; sigma[i] = bool IsPrime(ll n){
i+1; if(n <= 11) return hard_coding; 4.14 Power Tower
} if(n % 2 == 0 || ... 3 5 7 11) return false; bool PowOverflow(ll a, ll b, ll c){
for(auto j : prime){ for(int p : {comments}) if(!MillerRabin(n, p)) return false; __int128_t res = 1;
if(i*j >= sz) break; return true; bool flag = false;
sp[i*j] = j; } for(; b; b >>= 1, a = a * a){
if(i % j == 0){ ull Rho(ull n){ if(a >= c) flag = true, a %= c;
e[i*j] = e[i]+1; phi[i*j] = phi[i]*j; mu[i*j] = 0; ull x = 0, y = 0, t = 30, prd = 2, i = 1, q; if(b & 1){
tau[i*j] = tau[i]/e[i*j]*(e[i*j]+1); auto f = [&](ull x) { return MulMod(x, x, n) + i; }; res *= a; if(flag || res >= c) return true;
sigma[i*j] = sigma[i]*(j-1)/(pw(j, e[i*j])-1)*(pw(j, while(t++ % 40 || __gcd(prd, n) == 1){ }
e[i*j]+1)-1)/(j-1);//overflow if(x == y) x = ++i, y = f(x); }
break; if((q = MulMod(prd, max(x,y) - min(x,y), n))) prd = q; return false;
} x = f(x), y = f(f(y)); }
e[i*j] = 1; phi[i*j] = phi[i] * phi[j]; mu[i*j] = mu[i] * } return __gcd(prd, n); ll Recursion(int idx, ll mod, const vector<ll> &vec){
mu[j]; } if(mod == 1) return 1;
tau[i*j] = tau[i] * tau[j]; sigma[i*j] = sigma[i] * vector<ull> Factorize(ull n){ // sort? if(idx + 1 == vec.size()) return vec[idx];
sigma[j]; if(n == 1) return {}; if(IsPrime(n)) return {n}; ll nxt = Recursion(idx+1, phi[mod], vec);
} auto x = Rho(n); auto l=Factorize(x), r=Factorize(n/x); if(PowOverflow(vec[idx], nxt, mod)) return Pow(vec[idx], nxt,
} l.insert(l.end(), r.begin(), r.end()); return l; } mod) + mod; else return Pow(vec[idx], nxt, mod);
}
4.11 Xudyh Sieve 4.13 Primitive Root, Discrete Log/Sqrt ll PowerTower(const vector<ll> &vec, ll mod){ //
/* e(x) = [x==1], 1(x) = 1, id_k(x) = x^k ll PrimitiveRoot(ll p){ // order p-1 vec[0]^(vec[1]^(vec[2]^(...)))
mu: mobius function, id(x) = x vector<pair<ll,ll>> v = Factorize(p-1); if(vec.size() == 1) return vec[0] % mod;
Soongsil University – PS akgwi Page 18 of 25
else return Pow(vec[0], Recursion(1, phi[mod], vec), mod); } /* NTT : ang = pow(w, (mod-1)/n) % mod, inv_fft -> ang^{-1},
} if(r == -1) return false; root[i] = root[i-1] * ang
pivot(r, s); XOR Convolution: roots[*]=1, a[j+k] = u+v, a[j+k+i/2] = u-v
4.15 De Bruijn Sequence } OR Convolution: roots[*]=1, a[j+k+i/2] += inv_fft ? -u : u;
} AND Convolution: roots[*]=1, a[j+k ] += inv_fft ? -v : v; */
// Create cyclic string of length k^n that contains every length
// Returns -inf if no solution, {inf, a vector satisfying the for(int i=2; i<=N; i<<=1){ int step = N / i;
n string as substring. alphabet = [0, k - 1]
constraints} for(int j=0; j<N; j+=i) for(int k=0; k<i/2; k++){
int res[10000000], aux[10000000]; // >= k^n
// if there are abritrarily good solutions, or {maximum c^T cpx u = a[j+k], v = a[j+k+i/2] * root[step * k];
int de_bruijn(int k, int n) { // Returns size (k^n)
x, x} otherwise. a[j+k] = u+v; a[j+k+i/2] = u-v;
if(k == 1) { res[0] = 0; return 1; }
// O(n m (# of pivots)), O(2 ^ n) in general. } } // inv_fft: skip for AND/OR convolution.
for(int i = 0; i < k * n; i++) aux[i] = 0;
pair<T, vector<T>> solve(){ if(inv_fft) for(int i=0; i<N; i++) a[i] /= N;
int sz = 0;
int r = 0; }
function<void(int, int)> db = [&](int t, int p) {
for(int i=1; i<m; i++) if(mat[i][n+1] < mat[r][n+1]) r = i; vector<ll> Mul(const vector<ll> &_a, const vector<ll> &_b){
if(t > n) {
if(mat[r][n+1] < -eps){ vector<cpx> a(all(_a)), b(all(_b)); // (NTT) 2^19 700ms
if(n % p == 0) for(int i=1; i<=p; i++) res[sz++]=aux[i];
pivot(r, n); int N = 2; while(N < a.size() + b.size()) N <<= 1;
}
if(!simplex(2) || mat[m+1][n+1] < -eps) return {-inf, {}}; a.resize(N); b.resize(N); FFT(a); FFT(b);
else {
for(int i=0; i<m; i++) if(bb[i] == -1){ for(int i=0; i<N; i++) a[i] *= b[i]; // mod?
aux[t] = aux[t - p]; db(t + 1, p);
int s = 0; vector<ll> ret(N); FFT(a, 1); // NTT : just return a
for(int i=aux[t-p]+1; i<k; i++) aux[t]=i, db(t+1, t);
for(int j=1; j<n+1; j++) if(s == -1 || pair(mat[i][j], for(int i=0; i<N; i++) ret[i] = llround(a[i].real());
}
nn[j]) < pair(mat[i][s], nn[s])) s = j; while(ret.size() > 1 && ret.back() == 0) ret.pop_back();
}; db(1, 1); return sz;
pivot(i, s); return ret; }
}
} vector<ll> MulMod(const vector<ll> &a, const vector<ll> &b,
} const unsigned long long mod){ // (FFT) 2^19 1000ms
4.16 Simplex / LP Duality bool ok = simplex(1); int N = 2; while(N < a.size() + b.size()) N <<= 1;
// Solves the canonical form: maximize c^T x, subject to ax <= b vector<T> x(n); vector<cpx> v1(N), v2(N), r1(N), r2(N);
and x >= 0. for(int i=0; i<m; i++) if(bb[i] < n) x[bb[i]] = mat[i][n + for(int i=0;i<a.size();i++)v1[i] = cpx(a[i]>>15, a[i]&32767);
template<class T> // T must be of floating type 1]; for(int i=0;i<b.size();i++)v2[i] = cpx(b[i]>>15, b[i]&32767);
struct linear_programming_solver_simplex{ return {ok ? mat[m][n + 1] : inf, x}; FFT(v1); FFT(v2);
int m, n; vector<int> nn, bb; vector<vector<T>> mat; } for(int i=0; i<N; i++){ int j = i ? N-i : i;
static constexpr T eps = 1e-8, inf = 1/.0; }; cpx ans1 = (v1[i] + conj(v1[j])) * cpx(0.5, 0);
linear_programming_solver_simplex(const vector<vector<T>> &a, cpx ans2 = (v1[i] - conj(v1[j])) * cpx(0, -0.5);
const vector<T> &b, const vector<T> &c) : m(b.size()), Simplex Example cpx ans3 = (v2[i] + conj(v2[j])) * cpx(0.5, 0);
n(c.size()), nn(n+1), bb(m), mat(m+2, vector<T>(n+2)){ Maximize p = 6x + 14y + 13z cpx ans4 = (v2[i] - conj(v2[j])) * cpx(0, -0.5);
for(int i=0; i<m; i++) for(int j=0; j<n; j++) mat[i][j] = Constraints r1[i] = (ans1 * ans3) + (ans1 * ans4) * cpx(0, 1);
a[i][j]; - 0.5x + 2y + z ≤ 24 r2[i] = (ans2 * ans3) + (ans2 * ans4) * cpx(0, 1);
for(int i=0; i<m; i++) bb[i] = n + i, mat[i][n] = -1, - x + 2y + 4z ≤ 60 } vector<ll> ret(N); FFT(r1, true); FFT(r2, true);
mat[i][n + 1] = b[i]; Coding for(int i=0; i<N; i++){
for(int j=0; j<n; j++) nn[j] = j, mat[m][j] = -c[j]; 0.5 2 1 24 ll av = llround(r1[i].real()) % mod;
nn[n] = -1; mat[m + 1][n] = 1; - n = 2, m = 3, a = ,b = , c = [6, 14, 13]
1 2 4 60 ll bv = (llround(r1[i].imag()) + llround(r2[i].real()))%mod;
} LP Duality & Example ll cv = llround(r2[i].imag()) % mod;
void pivot(int r, int s){ tableu를 대각선으로 뒤집고 음수 부호를 붙인 답 =-(원문제의 답)
ret[i] = (av << 30) + (bv << 15) + cv;
T *a = mat[r].data(), inv = 1 / a[s]; 0.5 2 1 24 ret[i] %= mod; ret[i] += mod; ret[i] %= mod;
for(int i=0; i<m+2; i++) if(i != r && abs(mat[i][s]) > eps) - Primal : n = 2, m = 3, a = ,b = , c = [6, 14, 13] } while(ret.size() > 1 && ret.back() == 0) ret.pop_back();
1 2 4 60
{
−0.5 −1
−6
return ret; }
T *b = mat[i].data(), inv2 = b[s] * inv; - Dual : n = 3, m = 2, a = −2 −2 , b = −14 , c = [−24, −60] template<char op>vector<ll>FWHT_Conv(vector<ll> a,vector<ll> b){
for(int j=0; j<n+2; j++) b[j] -= a[j] * inv2; −1 −4 −13 int n = max({(int)a.size(), (int)b.size()-1, 1});//2^20 700ms
b[s] = a[s] * inv2; 공식 if(__builtin_popcount(n) != 1) n = 1 << (__lg(n) + 1);
} - Primal : maxx cT x, Constraints Ax ≤ b, x ≥ 0 a.resize(n); b.resize(n); FWHT<op>(a); FWHT<op>(b);
for(int j=0; j<n+2; j++) if(j != s) mat[r][j] *= inv; - Dual : miny bT y,Constraints AT y ≥ c, y ≥ 0 for(int i=0; i<n; i++) a[i] = a[i] * b[i] % M;
for(int i=0; i<m+2; i++) if(i != r) mat[i][s] *= -inv; FWHT<op>(a, true); return a;
mat[r][s] = inv; swap(bb[r], nn[s]); 4.17 Polynomial & Convolution } // subset: C[k] = sum_{i and j = 0, i or j = k} A[i] * B[j]
} // 998,244,353 = 119 23, w 3 | 2,281,701,377 = 17 27, w 3 vector<ll> SubsetConvolution(vector<ll> p,vector<ll> q){//Nlog2N
bool simplex(int phase){ // 167,772,161 = 10 25, w 3 | 2,483,027,969 = 37 26, w 3 int n = max({(int)p.size(), (int)q.size()-1, 1}), w=__lg(n);
for(auto x=m+phase-1; ; ){ // 469,762,049 = 26 26, w 3 | 2,013,265,921 = 15 27, w 31 if(__builtin_popcount(n) != 1) n = 1 << (w + 1); // 2^20 4s
int s = -1, r = -1; using real_t = double; using cpx = complex<real_t>; p.resize(n); q.resize(n); vector<ll> res(n); // SOS DP: 2.5s
for(auto j=0; j<n+1; j++) if(nn[j] != -phase) if(s == -1 void FFT(vector<cpx> &a, bool inv_fft=false){ vector<vector<ll>> a(w+1, vector<ll>(n)), b(a);
|| pair(mat[x][j], nn[j]) < pair(mat[x][s], nn[s])) s = j; int N = a.size(); vector<cpx> root(N/2);//root[0]=1 for(int i=0; i<n; i++) a[__builtin_popcount(i)][i] = p[i];
if(mat[x][s] >= -eps) return true; for(int i=1, j=0; i<N; i++){ int bit = N / 2; for(int i=0; i<n; i++) b[__builtin_popcount(i)][i] = q[i];
for(auto i=0; i<m; i++){ while(j >= bit) j -= bit, bit >>= 1; for(int bit=0; bit<=w; bit++) FWHT<'|'>(a[bit]),
if(mat[i][s] <= eps) continue; if(i < (j += bit)) swap(a[i], a[j]); FWHT<'|'>(b[bit]);
if(r == -1 || pair(mat[i][n + 1] / mat[i][s], bb[i]) < } long double ang = 2 * acosl(-1) / N * (inv_fft ? -1 : 1); for(int bit=0; bit<=w; bit++){
pair(mat[r][n + 1] / mat[r][s], bb[r])) r = i; for(int i=0; i<N/2; i++)root[i]=cpx(cosl(ang*i),sinl(ang*i));
Soongsil University – PS akgwi Page 19 of 25
vector<ll> c(n); // Warning : MOD return Trim(Integrate(Mul(Derivative(a), Inv(a,sz))), sz); } return a; }
for(int i=0; i<=bit; i++) for(int j=0; j<n; j++) c[j] += vector<ll> Exp(const vector<ll> &a, size_t sz){ // 5e5 5s vector<double> interpolate(vector<double> x, vector<double> y,
a[i][j] * b[bit-i][j] % M; vector<ll> res = {1}; if(a.empty()) return {1}; int n){ // n^2
for(auto &i : c) i %= M; assert(a.size() > 0 && a[0] == 0); vector<double> res(n), temp(n);
FWHT<'|'>(c, true); for(int i=1; i<sz; i<<=1){ for(int k=0; k<n-1; k++) for(int i=k+1; i<n; i++) y[i] = (y[i]
for(int i=0; i<n; i++) if(__builtin_popcount(i) == bit) auto t = Trim(a, i*2) - Log(res, i*2); - y[k]) / (x[i] - x[k]);
res[i] = c[i]; if(++t[0] == M) t[0] = 0; // t[0] += 1, mod double last = 0; temp[0] = 1;
} return res; } res = Trim(Mul(res, t), i*2); for(int k=0; k<n; k++){
vector<ll> Trim(vector<ll> a, size_t sz){ a.resize(min(a.size(), } return Trim(res, sz); } // need resize for(int i=0; i<n; i++) res[i] += y[k] * temp[i], swap(last,
sz)); return a; } vector<ll> Pow(const vector<ll> f, ll e, int sz){ // 5e5 8s temp[i]), temp[i] -= last * x[k];
vector<ll> Inv(const vector<ll> &a, size_t sz){ // 5e5 2s if(e == 0){ vector<ll> res(sz); res[0] = 1; return res; } }
vector<ll> q(1, Pow(a[0], M-2)); // 1/a[0], a[0] != 0 ll p = 0; while(p < f.size() && f[p] == 0 && p*e < sz) p++; return res; }//for numerical precision, x[k]=c*cos(k*pi/(n-1))
for(int i=1; i<sz; i<<=1){ // - : polynomial minus if(p == f.size() || p*e >= sz) return vector<ll>(sz, 0); vector<ll> Interpolation_0_to_n(vector<ll> y){ // n^2
auto p = vector<ll>{2} - Mul(q, Trim(a, i*2)); vector<ll> a(f.begin()+p, f.end()); ll k = a[0]; // not f[0] int n = y.size();
q = Trim(Mul(p, q), i*2); for(auto &i : a) i = mul(i, Pow(k, M-2)); vector<ll> res(n), tmp(n), x; // x[i] = i / (i+1)
} return Trim(q, sz); } a = Log(a, sz); for(auto &i : a) i = mul(i, e%M); for(int i=0; i<n; i++) x.push_back(Pow(i+1, M-2));
vector<ll> Div(const vector<ll> &a, const vector<ll> &b){ a = Exp(a, sz); for(auto &i : a) i = mul(i, Pow(k, e)); for(int k=0; k+1<n; k++) for(int i=k+1; i<n; i++)
if(a.size() < b.size()) return {}; // 5e5 4s vector<ll> res(p*e); res.insert(res.end(), a.begin(), y[i] = (y[i] - y[k] + M) * x[i-k-1] % M;
size_t sz = a.size() - b.size() + 1; auto ra = a, rb = b; a.end()); ll lst = 0; tmp[0] = 1;
reverse(ra.begin(), ra.end()); ra = Trim(ra, sz); res.resize(sz); return res; } for(int k=0; k<n; k++) for(int i=0; i<n; i++) {
reverse(rb.begin(), rb.end()); rb = Inv(Trim(rb,sz), sz); vector<ll> SqrtImpl(vector<ll> a){ res[i] = (res[i] + y[k] * tmp[i]) % M;
auto res = Trim(Mul(ra, rb), sz); res.resize(sz); if (a.empty()) return {0}; int inv2=(M+1)/2; swap(lst, tmp[i]);
reverse(res.begin(), res.end()); int z = DiscreteSqrt(a[0], M), n = a.size(); tmp[i] = (tmp[i] - lst * k) % M;
while(!res.empty() && !res.back()) res.pop_back(); if (z == -1) return {-1}; vector<ll> q(1, z); if(tmp[i] < 0) tmp[i] += M;
return res; } for(int m=1; m<n; m<<=1){ } return res; }
vector<ll> Mod(const vector<ll> &a, const vector<ll> &b){ return if(n < m*2) a.resize(m*2);;; q.resize(m*2);
a - Mul(b, Div(a, b)); } auto f2 = Mul(q, q); f2.resize(m*2);
ll Evaluate(const vector<ll> &a, ll x){ ll res = 0; for(int i=0; i<m*2; i++) f2[i] = sub(f2[i], a[i]);
4.18 Matroid Intersection
for(int i=(int)a.size()-1; i>=0; i--) res = (res*x+a[i]) % M; f2 = Mul(f2, Inv(q, q.size())); f2.resize(m*2); struct Matroid{
return res >= 0 ? res : res + M; } for(int i=0; i<m*2; i++) q[i] = sub(q[i], mul(f2[i], inv2)); virtual bool check(int i) = 0; // O(R^2N), O(R^2N)
vector<ll> Derivative(const vector<ll> &a){ } q.resize(n); return q; } virtual void insert(int i) = 0; // O(R^3), O(R^2N)
if(a.size() <= 1) return {}; vector<ll> res(a.size()-1); vector <ll> Sqrt(vector <ll> a){ // nlgn, fail -> -1, 5e5 5.5s virtual void clear() = 0; // O(R^2), O(RN)
for(int i=0; i+1<a.size(); i++) res[i] = (i+1) * a[i+1] % M; int n = a.size(), m = 0; while(m < n && a[m] == 0) m++; };
return res; } if(m == n) return vector<ll>(n); if(m & 1) return {-1}; template<typename cost_t>
vector<ll> Integrate(const vector<ll> &a){ auto s = SqrtImpl(vector<ll>(a.begin()+m, a.end())); vector<cost_t> MI(const vector<cost_t> &cost, Matroid *m1,
int n = a.size(); vector<ll> res(n+1); if(s[0] == -1) return {-1}; vector<ll> res(n); Matroid *m2){
for(int i=0; i<n; i++) res[i+1] = a[i] * Pow(i+1, M-2) % M; for(int i=0; i<s.size(); i++) res[i+m/2] = s[i]; int n = cost.size();
return res; } return res; } vector<pair<cost_t, int>> dist(n+1);
vector<ll> MultipointEvaluation(vector<ll> a, vector<ll> x){ vector<ll> TaylorShift(vector<ll> a, ll c){//f(x+c), 2^19 700ms vector<vector<pair<int, cost_t>>> adj(n+1);
if(x.empty()) return {}; int n = x.size(); // 2^17 7s int n = a.size(); // fac[i] = i!, ifc[i] = inv(i!) vector<int> pv(n+1), inq(n+1), flag(n); deque<int> dq;
vector<vector<ll>> up(n*2), dw(n*2); for(int i=0; i<n; i++) a[i] = mul(a[i], fac[i]); auto augment = [&]() -> bool {
for(int i=0; i<n; i++) up[i+n] = {x[i]?M-x[i]:0, 1}; reverse(all(a)); vector<ll> b(n); ll w = 1; fill(dist.begin(), dist.end(),
for(int i=n-1; i; i--) up[i] = Mul(up[i*2], up[i*2+1]); for(int i=0; i<n; i++) b[i] = mul(ifc[i], w), w = mul(w, c); pair(numeric_limits<cost_t>::max()/2, 0));
dw[1] = Mod(a, up[1]); a = Mul(a, b); a.resize(n); reverse(all(a)); fill(adj.begin(), adj.end(), vector<pair<int, cost_t>>());
for(int i=2; i<n*2; i++) dw[i] = Mod(dw[i/2], up[i]); for(int i=0; i<n; i++) a[i] = mul(a[i], ifc[i]); fill(pv.begin(),pv.end(),-1); fill(inq.begin(),inq.end(),0);
vector<ll> y(n); for(int i=0; i<n; i++) y[i] = dw[i+n][0]; return a; } dq.clear(); m1->clear(); m2->clear();
return y; } vector<ll> SamplingShift(vector<ll> a, ll c, int m){ // 2^19 ~2s for(int i=0;i<n;i++)if(flag[i])m1->insert(i),m2->insert(i);
vector<ll> Interpolation(vector<ll> x, vector<ll> y){//2^17 10s // given f(0), f(1), ..., f(n - 1), warning: fac size for(int i=0; i<n; i++){
int n = x.size(); vector<vector<ll>> up(n*2), dw(n*2); // return f(c), f(c + 1), ..., f(c + m - 1) if(flag[i]) continue;
for(int i=0; i<n; i++) up[i+n] = {x[i]?M-x[i]:0, 1}; int n = a.size(); vector<ll> b(ifc.begin(), ifc.begin()+n); if(m1->check(i))
for(int i=n-1; i; i--) up[i] = Mul(up[i*2], up[i*2+1]); for(int i=0; i<n; i++) a[i] = mul(a[i],ifc[i]); dist[pv[i]=i] = {cost[i], 0}, dq.push_back(i), inq[i]=1;
vector<ll> a = MultipointEvaluation(Derivative(up[1]), x); for(int i=1; i<n; i+=2) b[i] = sub(0, b[i]); if(m2->check(i)) adj[i].emplace_back(n, 0);
for(int i=0; i<n; i++) a[i] = y[i] * Pow(a[i], M-2) % M; a = Mul(a, b); a.resize(n); ll w = 1; }
for(int i=0; i<n; i++) dw[i+n] = {a[i]}; for(int i=0; i<n; i++) a[i] = mul(a[i], fac[i]);; for(int i=0; i<n; i++){
for(int i=n-1; i; i--){ reverse(all(a)); if(!flag[i]) continue; m1->clear(); m2->clear();
auto l = Mul(dw[i*2],up[i*2+1]), r = Mul(dw[i*2+1],up[i*2]); for(int i=0; i<n; w=mul(w, sub(c,i++))) b[i] = mul(ifc[i],w); for(int j=0; j<n; j++) if(i != j && flag[j])
dw[i].resize(l.size()); a = Mul(a, b); a.resize(n); reverse(all(a)); //warning: N->M m1->insert(j), m2->insert(j);
for(int j=0; j<l.size(); j++) dw[i][j] = (l[j] + r[j]) % M; for(int i=0; i<n; i++) a[i] = mul(a[i], ifc[i]);; a.resize(m); for(int j=0; j<n; j++){
} return dw[1]; } b = vector<ll>(ifc.begin(), ifc.begin()+m); if(flag[j]) continue;
vector<ll> Log(const vector<ll> &a, size_t sz){ // 5e5 3.5s a = Mul(a, b); a.resize(m); if(m1->check(j)) adj[i].emplace_back(j, cost[j]);
assert(a.size() > 0 && a[0] == 1); // int f'(x)/f(x), resize! for(int i=0; i<m; i++) a[i] = mul(a[i], fac[i]); if(m2->check(j)) adj[j].emplace_back(i, -cost[i]);
}
Soongsil University – PS akgwi Page 20 of 25
} int n = s.size(); vector<int> p(n); auto counting_sort = [&](){
while(dq.size()){ vector<pair<int,int>> range, maximal; fill(cnt.begin(), cnt.end(), 0);
int v = dq.front(); dq.pop_front(); inq[v] = 0; auto make = [&](int l, int r) { return make_pair(l/2, for(int i=0; i<n; i++) cnt[pos[i]]++;
for(const auto &[i,w] : adj[v]){ (r-1)/2); }; partial_sum(cnt.begin(), cnt.end(), cnt.begin());
pair<cost_t, int> nxt{dist[v].ff+w, dist[v].ss+1}; for(int i=0, k=-1, r=-1; i<n; i++){ for(int i=n-1; i>=0; i--) sa[--cnt[pos[tmp[i]]]] = tmp[i];
if(nxt < dist[i]){ if(i <= r) p[i] = min(r-i, p[2*k-i]); };
dist[i] = nxt; pv[i] = v; while(i-p[i]-1 >= 0 && i+p[i]+1<n && s[i-p[i]-1] == for(int i=0; i<n; i++) sa[i] = i, pos[i] = s[i], tmp[i] = i;
if(!inq[i]) dq.push_back(i), inq[i] = 1; s[i+p[i]+1]){ counting_sort();
} /* if */ } /* for [i,w] */ } /* while */ p[i]++; range.push_back(make(i-p[i], i+p[i])); for(int k=1; ; k<<=1){ int p = 0;
if(pv[n] == -1) return false; } if(i+p[i] > r) r = i+p[i], k = i; for(int i=n-k; i<n; i++) tmp[p++] = i;
for(int i=pv[n]; ; i=pv[i]){ if(p[i] != 0) maximal.push_back(make(i-p[i], i+p[i])); for(int i=0; i<n; i++) if(sa[i] >= k) tmp[p++] = sa[i] - k;
flag[i] ^= 1; if(i == pv[i]) break; } // compress(range), range can contains O(1) dup. substr... counting_sort(); tmp[sa[0]] = 0;
} return true; return p; } // range: distinct palindrome(<= n) for(int i=1; i<n; i++){
}; vector<cost_t> res; //z[i]=match length of s[0,n-1] and s[i,n-1] tmp[sa[i]] = tmp[sa[i-1]];
while(augment()){ vector<int> Z(const string &s){ if(sa[i-1]+k < n && sa[i]+k < n && pos[sa[i-1]] ==
cost_t now = cost_t(0); int n = s.size(); vector<int> z(n); z[0] = n; pos[sa[i]] && pos[sa[i-1]+k] == pos[sa[i]+k]) continue;
for(int i=0; i<n; i++) if(flag[i]) now += cost[i]; for(int i=1, l=0, r=0; i<n; i++){ tmp[sa[i]] += 1;
res.push_back(now); if(i < r) z[i] = min(r-i-1, z[i-l]); }
} return res; while(i+z[i] < n && s[i+z[i]] == s[z[i]]) z[i]++; swap(pos, tmp); if(pos[sa.back()] + 1 == n) break;
} if(i+z[i] > r) r = i+z[i], l = i; }
} return z; for(int i=0, j=0; i<n; i++, j=max(j-1,0)){
5 String } if(pos[i] == 0) continue;
while(sa[pos[i]-1]+j < n && sa[pos[i]]+j < n &&
5.1 KMP, Hash, Manacher, Z
5.2 Aho-Corasick s[sa[pos[i]-1]+j] == s[sa[pos[i]]+j]) j++;
vector<int> getFail(const container &pat){ lcp[pos[i]] = j;
vector<int> fail(pat.size()); struct Node{
int g[26], fail, out; } return {sa, lcp};
//match: pat[0..j] and pat[j-i..i] is equivalent }
//ins/del: manipulate corresponding range to pattern starts at 0 Node() { memset(g, 0, sizeof g); fail = out = 0; }
}; auto [SA,LCP] = SuffixArray(S); RMQ<int> rmq(LCP);
// (insert/delete pat[i], manage pat[j-i..i]) vector<int> Pos(N); for(int i=0; i<N; i++) Pos[SA[i]] = i;
function<bool(int, int)> match = [&](int i, int j){ }; vector<Node> T(2); int aut[100101][26];
void Insert(int n, int i, const string &s){ auto get_lcp = [&](int a, int b){
function<void(int)> ins = [&](int i){ }; if(Pos[a] > Pos[b]) swap(a, b);
function<void(int)> del = [&](int i){ }; if(i == s.size()){ T[n].out++; return; }
int c = s[i] - 'a'; return a == b ? (int)S.size() - a : rmq.query(Pos[a]+1,
for(int i=1, j=0; i<pat.size(); i++){ Pos[b]);
while(j && !match(i, j)){ if(T[n].g[c] == 0) T[n].g[c] = T.size(), T.emplace_back();
Insert(T[n].g[c], i+1, s); };
for(int s=i-j; s<i-fail[j-1]; s++) del(s); vector<pair<int,int>> can; // common substring {start, lcp}
j = fail[j-1]; } }
int go(int n, int i){ // DO NOT USE `aut` DIRECTLY vector<tuple<int,int,int>> valid; // valid substring [string,
if(match(i, j)) ins(i), fail[i] = ++j; end_l~end_r]
} return fail; int &res = aut[n][i]; if(res) return res;
if(n != 1 && T[n].g[i] == 0) res = go(T[n].fail, i); for(int i=1; i<N; i++){
} if(SA[i] < X && SA[i-1] > X) can.emplace_back(SA[i], LCP[i]);
vector<int> doKMP(const container &str, const container &pat){ else if(T[n].g[i] != 0) res = T[n].g[i]; else res = 1;
return res; if(i+1 < N && SA[i] < X && SA[i+1] > X)
vector<int> ret, fail = getFail(pat); can.emplace_back(SA[i], LCP[i+1]);
//match: pat[0..j] and str[j-i..i] is equivalent }
void Build(){ }
//ins/del: manipulate corresponding range to pattern starts at 0 for(int i=0; i<can.size(); i++){
// (insert/delete str[i], manage str[j-i..i]) queue<int> q; q.push(1); T[1].fail = 1;
while(!q.empty()){ int skip = i > 0 ? min({can[i-1].second, can[i].second,
function<bool(int, int)> match = [&](int i, int j){ }; get_lcp(can[i-1].first, can[i].first)}) : 0;
function<void(int)> ins = [&](int i){ }; int n = q.front(); q.pop();
for(int i=0; i<26; i++){ valid.emplace_back(can[i].first, can[i].first + skip,
function<void(int)> del = [&](int i){ }; can[i].first + can[i].second - 1);
for(int i=0, j=0, s; i<str.size(); i++){ int next = T[n].g[i]; if(next == 0) continue;
if(n == 1)T[next].fail=1;else T[next].fail=go(T[n].fail,i); }
while(j && !match(i, j)){
for(int s=i-j; s<i-fail[j-1]; s++) del(s); q.push(next); T[next].out += T[T[next].fail].out;
j = fail[j-1]; } } /* for i */ } /* while q */ } /* build */
if(match(i, j)){ bool Find(const string &s){ 5.4 O(N log N ) Tandem Repeats
if(j+1 == pat.size()){ int n = 1, ok = 0; // return O(n log n) tuple {l, r, p} that
ret.push_back(i-j); for(s=i-j;s<i-fail[j]+1; s++)del(s); for(int i=0; i<s.size(); i++){ // [i, i+p) = [i+p, i+2p) for all l <= i < r
j = fail[j]; n = go(n, s[i] - 'a'); if(T[n].out != 0) ok = 1; vector<tuple<int,int,int>> TandemRepeat(const string &s){
} else ++j; ins(i); } return ok; int n = s.size(); vector<tuple<int,int,int>> res;
} } return ret; } string t = s; reverse(t.begin(), t.end());
} // WARNING: add empty suffix!!
// # a # b # a # a # b # a # 5.3 O(N log N ) SA + LCP // sa.insert(sa.begin(), n) before calculate lcp/pos
// 0 1 0 3 0 1 6 1 0 3 0 1 0 pair<vector<int>, vector<int>> SuffixArray(const string &s){ auto [sa_s,lcp_s,pos_s] = SuffixArray(s);
vector<int> Manacher(const string &inp){ int n = s.size(), m = max(n, 256); auto [sa_t,lcp_t,pos_t] = SuffixArray(t);
string s = "#"; for(auto c : inp) s += c, s += "#"; vector<int> sa(n), lcp(n), pos(n), tmp(n), cnt(m); RMQ<int> rmq_s(lcp_s), rmq_t(lcp_t);
Soongsil University – PS akgwi Page 21 of 25
auto get = [n](const vector<int> &pos, const RMQ<int> &rmq, 5.6 Bitset LCS 6 Misc
int a, int b){ #include <x86intrin.h> 6.1 CMakeLists.txt
if(pos[a] > pos[b]) swap(a, b); template<size_t _Nw> void _M_do_sub(_Base_bitset<_Nw> &A, const
return a == b ? n - a : rmq.query(pos[a] + 1, pos[b]); set(CMAKE_CXX_STANDARD 17)
_Base_bitset<_Nw> &B){
}; set(CMAKE_CXX_FLAGS "-DLOCAL -lm -g -Wl,--stack,268435456")
for(int i=0, c=0; i<_Nw; i++) c = _subborrow_u64(c, A._M_w[i],
for(int p=1; p*2<=n; p++){ add_compile_options(-Wall -Wextra -Winvalid-pch -Wfloat-equal
B._M_w[i], (ull*)&A._M_w[i]);
for(int i=0, j=-1; i+p<=n; i+=p){ -Wno-sign-compare -Wno-misleading-indentation -Wno-parentheses)
}
int l = i - get(pos_t, rmq_t, n-i-p, n-i); # add_compile_options(-O3 -mavx -mavx2 -mfma)
void _M_do_sub(_Base_bitset<1> &A, const _Base_bitset<1> &B){
int r = i - p + get(pos_s, rmq_s, i, i+p); A._M_w -= B._M_w; }
if(l <= r && l != j) res.emplace_back(j=l, r+1, p); template<size_t _Nb> bitset<_Nb>& operator-=(bitset<_Nb> &A, 6.2 Calendar
}} return res; const bitset<_Nb> &B){ int f(int y,int m,int d){// 0: Sat, 1: Sun, ...
} // Check p = 0, time complexity O(n log n) _M_do_sub(A, B); return A; if (m<=2) y--, m+=12; int c=y/100; y%=100;
5.5 Suffix Automaton } int w=((c>>2)-(c<<1)+y+(y>>2)+(13*(m+1)/5)+d-1)%7;
template<size_t _Nb> inline bitset<_Nb> operator-(const if (w<0) w+=7; return w; }
template<typename T, size_t S, T init_val>
bitset<_Nb> &A, const bitset<_Nb> &B){
struct initialized_array : public array<T, S> {
bitset<_Nb> C(A); return C -= B;
initialized_array(){ this->fill(init_val); } 6.3 Ternary Search
}
};
char s[50050], t[50050]; while(s + 3 <= e){
template<class Char_Type, class Adjacency_Type>
int lcs(){ // O(NM/64) T l = (s + s + e) / 3, r = (s + e + e) / 3;
struct suffix_automaton{
bitset<50050> dp, ch[26]; if(Check(l) > Check(r)) s = l; else e = r;
// Begin States
int n = strlen(s), m = strlen(t); }// get minimum / when multiple answer, find minimum `s`
// len: length of the longest substring in the class
for(int i=0; i<m; i++) ch[t[i]-'A'].set(i); T mn = INF, idx = s;
// link: suffix link
for(int i=0; i<n; i++){ auto x = dp | ch[s[i]-'A']; dp = dp - for(T i=s; i<=e; i++) if(T now = Check(i); now < mn) mn = now,
// firstpos: minimum value in the set endpos
(dp ^ x) & x; } idx = i;
vector<int> len{0}, link{-1}, firstpos{-1}, is_clone{false};
return dp.count();
vector<Adjacency_Type> next{{}};
}
ll ans{0LL}; // 서로 다른 부분 문자열 개수 6.4 Add/Mul Update, Range Sum Query
// End States
void set_link(int v, int lnk){
5.7 Lyndon Factorization, Minimum Rotation struct Lz{
// link[i]: length of smallest suffix of s[0..i-1] ll a, b; // constructor, clear(a = 1, b = 0)
if(link[v] != -1) ans -= len[v] - len[link[v]]; Lz& operator+=(const Lz &t); // a *= t.a, b = t.a * b + t.b
link[v] = lnk; // factorization result: s[res[i]..res[i+1]-1]
vector<int> Lyndon(const string &s){ };
if(link[v] != -1) ans += len[v] - len[link[v]]; struct Ty{
} int n = s.size(); vector<int> link(n);
for(int i=0; i<n; ){ ll cnt, sum; // constructor cnt=1, sum=0
int new_state(int l, int sl, int fp, bool c, const Ty& operator += (const Ty &t); // cnt += t.cnt, sum += t.sum
Adjacency_Type &adj){ int j=i+1, k=i; link[i] = 1;
for(; j<n && s[k]<=s[j]; j++){ Ty& operator += (const Lz &t); // sum= t.a * sum + cnt * t.b
int now = len.size(); len.push_back(l); link.push_back(-1); };
set_link(now, sl); firstpos.push_back(fp); if(s[j] == s[k]) link[j] = link[k], k++;
is_clone.push_back(c); next.push_back(adj); return now; else link[j] = j - i + 1, k = i;
} int last = 0; } for(; i<=k; i+=j-k); 6.5 O(N × max W ) Subset Sum (Fast Knapsack)
void extend(const vector<Char_Type> &s){ } vector<int> res;
// O(N*maxW), maximize sumW <= t
last = 0; for(auto c: s) extend(c); } for(int i=n-1; i>=0; i-=link[i]) res.push_back(i-link[i]+1);
int Knapsack(vector<int> w, int t){
void extend(Char_Type c){ reverse(res.begin(), res.end()); return res;
int a = 0, b = 0, x;
int cur = new_state(len[last] + 1, -1, len[last], false, }
while(b < w.size() && a + w[b] <= t) a += w[b++];
{}), p = last; // rotate(v.begin(), v.begin()+min_rotation(v), v.end());
if(b == w.size()) return a;
while(~p && !next[p][c]) next[p][c] = cur, p = link[p]; template<typename T> int min_rotation(T s){ // O(N)
int m = *max_element(w.begin(), w.end());
if(!~p) set_link(cur, 0); int a = 0, N = s.size();
vector<int> u, v(2*m, -1); v[a+m-t] = b;
else{ for(int i=0; i<N; i++) s.push_back(s[i]);
for(int i=b; (u=v,i<w.size()); i++){
int q = next[p][c]; for(int b=0; b<N; b++) for(int k=0; k<N; k++){
for(x=0; x<m; x++) v[x+w[i]] = max(v[x+w[i]], u[x]);
if(len[p] + 1 == len[q]) set_link(cur, q); if(a+k == b || s[a+k] < s[b+k]){ b += max(0, k-1); break; }
for(x=2*m; --x>m; ) for(int j=max(0,u[x]); j<v[x]; j++)
else{ if(s[a+k] > s[b+k]){ a = b; break; }
v[x-w[j]] = max(v[x-w[j]], j);
int clone = new_state(len[p] + 1, link[q], firstpos[q], }
} for(a=t; v[a+m-t]<0; a--);;; return a;
true, next[q]); return a;
}
while(~p && next[p][c] == q) next[p][c] = clone, p = }
link[p];
set_link(cur, clone); set_link(q, clone); 5.8 All LCS 6.6 Monotone Queue Optimization
} void AllLCS(const string &s, const string &t){ template<class T, bool GET_MAX = false> // D[i] = func_{0 <= j <
} vector<int> h(t.size()); iota(h.begin(), h.end(), 0); i} D[j] + cost(j, i)
last = cur; for(int i=0, v=-1; i<s.size(); i++, v=-1){ pair<vector<T>, vector<int>> monotone_queue_dp(int n, const
} for(int r=0; r<t.size(); r++){ vector<T> &init, auto cost){
int size() const { return (int)len.size(); } // # of states if(s[i] == t[r] || h[r] < v) swap(h[r], v); assert((int)init.size() == n + 1); // cost function -> auto,
}; suffix_automaton<int, initialized_array<int,26,0>> T; //LCS(s[0..i],t[l..r]) = r-l+1 - sum([h[x] >= l] | x <= r) do not use std::function
// for(auto c : s) if((x=T.next[x][c]) == 0) return false; } /*for r*/ } /* for i */ } /* end*/ vector<T> dp = init; vector<int> prv(n+1);
Soongsil University – PS akgwi Page 22 of 25
auto compare = [](T a, T b){ return GET_MAX ? a < b : a > b; //x[i]>ans[i+1]일때: ans[i]=ans[i+1] 왜냐하면 f(i,a)는 a<x[i]에서 for(int i=0; i<x; i+=2){
}; 감소함수이므로 가능한 최대로 오른쪽으로 붙은 ans[i+1]이 최적. int to = i+1 < x ? ans[rs[i+1]] : cs[y-1]; T mx;
auto cross = [&](int i, int j){ //스텝i에서 add_bias(k,0)한다면 간격제한k가 있는것이므로 while(true){
int l = j, r = n + 1; ans[i]=min(ans[i+1]-k,x[i])으로 수정. T w = f(rs[i], cs[k]);
while(l < r){ //LR Hull 역추적은 케이스나눠서 위 방법을 확장하면 될듯 if(ans[rs[i]] == -1 || w > mx) ans[rs[i]]=cs[k], mx=w;
int m = (l + r + 1) / 2; if(cs[k] == to) break; k++;
if(compare(dp[i] + cost(i, m), dp[j] + cost(j, m))) r = m 6.8 Aliens Trick }}
- 1; else l = m; // pair<T, vector<int>> f(T c): return opt_val, prv };
} return l; }; // cost function must be multiplied by 2 int *rs = new int[n*2]; iota(rs,rs+n,0);
deque<int> q{0}; template<class T, bool GET_MAX = false> int *cs = new int[max(m, n*2)]; iota(cs,cs+m,0);
for(int i=1; i<=n; i++){ pair<T,vector<int>> AliensTrick(int n,int k,auto f,T lo,T hi){ rec(rec,rs,n,cs,m);delete[]rs;delete[]cs;return ans;
while(q.size() > 1 && compare(dp[q[0]] + cost(q[0], i), T l = lo, r = hi; while(l < r) { }
dp[q[1]] + cost(q[1], i))) q.pop_front(); T m = (l + r + (GET_MAX?1:0)) >> 1; // A: convex, B: arbitrary, O((N+M) log M)
dp[i] = dp[q[0]] + cost(q[0], i); prv[i] = q[0]; vector<int> prv = f(m*2+(GET_MAX?-1:+1)).second; int N, M, A[1<<19], B[1<<19], C[1<<20];
while(q.size() > 1 && cross(q[q.size()-2], q.back()) >= int cnt = 0; for(int i=n; i; i=prv[i]) cnt++; void DnC(int s, int e, int l, int r){
cross(q.back(), i)) q.pop_back(); if(cnt <= k) (GET_MAX?l:r) = m; if(s > e) return;
q.push_back(i); else (GET_MAX?r:l) = m + (GET_MAX?-1:+1); int m = (s+e)/2, opt = -1, &mn = C[m];
} /*for end*/ return {dp, prv}; } } for(int i=l; i<=min(m,r); i++){
T opt_value = f(l*2).first / 2 - k*l; if(m - i >= N) continue;
6.7 Slope Trick if(opt == -1 || A[m-i] + B[i] < mn) mn=A[m-i]+B[i], opt=i;
//NOTE: f(x)=min{f(x+i),i<a}+|x-k|+m -> pf(k)sf(k)ab(-a,m) vector<int> prv1 = f(l*2+(GET_MAX?1:-1)).second, p1{n};
vector<int> prv2 = f(l*2-(GET_MAX?1:-1)).second, p2{n}; } DnC(s, m-1, l, opt); DnC(m+1, e, opt, r);
//NOTE: sf_inc에 답구하는게 들어있어서, 반드시 한 연산에 대해 } // or...
pf_dec->sf_inc순서로 호출 for(int i=n; i; i=prv1[i]) p1.push_back(prv1[i]);
for(int i=n; i; i=prv2[i]) p2.push_back(prv2[i]); int f(int r, int c){//O(N+M) but not fast
struct LeftHull{ if(0 <= r-c && r-c < N) return -(A[r-c] + B[c]);
void pf_dec(int x){ pq.empl(x-bias); }//x이하의 기울기들 -1 reverse(p1.begin(),p1.end());reverse(p2.begin(),p2.end());
assert(p2.size() <= k+1 && k+1 <=p1.size()); else return -21e8 - (r - c); // min
int sf_inc(int x){//x이상의 기울기들 +1, pop된 원소 반환(Right } SMAWK(f, N+M-1, M); // DnC opt 163ms SMAWK 179ms N,M=2^19
Hull관리에 사용됨) if(p1.size() == k+1) return {opt_value, p1};
if(pq.empty() or argmin()<=x) return x; ans += argmin()-x;// if(p2.size() == k+1) return {opt_value, p2};
이 경우 최솟값이 증가함 for(int i=1, j=1; i<p1.size(); i++){ 6.10 Money for Nothing (WF17D)
pq.empl(x-bias);/*x 이하 -1*/int r=argmin(); pq.pop();/*전체 while(j < p2.size() && p2[j] < p1[i-1]) j++;
ll MoneyForNothing(vector<Point> lo, vector<Point> hi){
+1*/ if(p1[i] <= p2[j] && i - j == k+1 - (int)p2.size()){
sort(lo.begin(), lo.end()); sort(hi.rbegin(), hi.rend());//rev
return r; vector<int> res;
vector<Point> a, b; ll res = 0;
} res.insert(res.end(), p1.begin(), p1.begin()+i);
for(auto p:lo)if(a.empty() || a.back().y > p.y)a.push_back(p);
void add_bias(int x,int y){ bias+=x; ans+=y; } int minval(){ res.insert(res.end(), p2.begin()+j, p2.end());
for(auto p:hi)if(b.empty() || b.back().y < p.y)b.push_back(p);
return ans; } //x축 평행이동, 최소값 return {opt_value, res};
reverse(b.begin(),b.end()); if(a.empty()||b.empty()) return 0;
int argmin(){return pq.empty()?-inf<int>():pq.top()+bias;}// } /* if */ } /* for */ assert(false);
queue<tuple<int,int,int,int>> que;
최소값 x좌표 }
que.emplace(0, (int)a.size()-1, 0, (int)b.size()-1);
void operator+=(LeftHull& a){ ans+=a.ans; while(sz(a.pq)) while(!que.empty()){
pf_dec(a.argmin()), a.pq.pop(); } 6.9 SWAMK, Min Plus Convolution auto [s,e,l,r] = que.front(); que.pop();
int size()const{return sz(pq);} PQMax<int> pq; int ans=0, // find the indices of row maxima, smallest index when tie int m = (s + e) / 2, pos = l; ll mx = -4e18;
bias=0; template <class F, class T=long long> for(int i=l; i<=r; i++){
}; vector<int> SMAWK(F f, int n, int m){ ll dx = b[i].x - a[m].x, dy = b[i].y - a[m].y;
//NOTE: f(x)=min{f(x+i),a<i<b}+|x-k|+m -> pf(k)sf(k)ab(-a,b,m) vector<int> ans(n, -1); ll now = (dx < 0 && dy < 0) ? 0 : dx * dy;
struct SlopeTrick{ auto rec = [&](auto self, int*const rs, int x, int*const cs, if(now > mx) mx = now, pos = i;
void pf_dec(int x){l.pf_dec(-r.sf_inc(-x));} int y){ } res = max(res, mx);
void sf_inc(int x){r.pf_dec(-l.sf_inc(x));} const int t = 8; if(s < m) que.emplace(s, m-1, l, pos);
void add_bias(int lx,int rx,int if(x <= t || y <= t){ if(m < e) que.emplace(m+1, e, pos, r);
y){l.add_bias(lx,0),r.add_bias(-rx,0),ans+=y;} for(int i=0; i<x; i++){ int r = rs[i]; T mx; } return res;
int minval(){return ans+l.minval()+r.minval();} for(int j=0; j<y; j++){ }
pint argmin(){return {l.argmin(),-r.argmin()};} int c = cs[j]; T w = f(r, c);
void operator+=(SlopeTrick& a){ if(ans[r] == -1 || w > mx) ans[r] = c, mx = w;
while(sz(a.l.pq)) pf_dec(a.l.argmin()),a.l.pq.pop(); }} /* for j i */ return; } /* base case */ 6.11 Exchange Argument on Tree (WF16L,CERC13E)
l.ans+=a.l.ans; if(x < y){ int s = 0; struct Info{ // down a -> up b, a b >= 0
while(sz(a.r.pq)) sf_inc(-a.r.argmin()),a.r.pq.pop(); for(int i=0; i<y; i++){ int c = cs[i]; ll a, b, idx; Info() : Info(0, 0, 0) {}
r.ans+=a.r.ans; ans+=a.ans; while(s && f(rs[s-1], cs[s-1]) < f(rs[s-1], c)) s--; Info(ll a, ll b, ll idx) : a(a), b(b), idx(idx) {}
} LeftHull l,r; int ans=0; if(s < x) cs[s++] = c; bool operator < (const Info &t) const {
int size()const{return l.size()+r.size();} } y = s; ll le = b - a, ri = t.b - t.a;
}; } if(le >= 0 && ri < 0) return false;
//LeftHull 역추적 방법: 스텝i의 argmin값을 am(i)라고 하자. 스텝n부터 int z=0, k=0, *a=rs+x, *b=cs+y; if(le < 0 && ri >= 0) return true;
스텝1까지 ans[i]=min(ans[i+1],am(i))하면 된다. 아래는 증명..은 for(int i=1; i<x; i+=2) a[z++] = rs[i]; if(le < 0 && b != t.b) return b < t.b;
아니고 간략한 이유 for(int i=0; i<y; i++) b[i] = cs[i]; if(le >= 0 && a != t.a) return a > t.a;
//am(i)<=ans[i+1]일때: ans[i]=am(i) self(self, a, z, b, y); return idx < t.idx;
Soongsil University – PS akgwi Page 23 of 25
} return r; } 11 97772875200 4032 6 3 2 2 1 1 1 1
Info& operator += (const Info &v){ bitset<17> bs; bs[1] = bs[7] = 1; assert(bs._Find_first() == 1); 12 963761198400 6720 6 4 2 1 1 1 1 1 1
ll aa = min(-a, -a+b-v.a), bb = -a+b-v.a+v.b; assert(bs._Find_next(0) == 1 && bs._Find_next(1) == 7); 13 9316358251200 10752 6 3 2 1 1 1 1 1 1 1
a = -aa; b = bb - aa; return *this; assert(bs._Find_next(3) == 7 && bs._Find_next(7) == 17); 14 97821761637600 17280 5 4 2 2 1 1 1 1 1 1
} cout << bs._Find_next(7) << "\n"; 15 866421317361600 26880 6 4 2 1 1 1 1 1 1 1 1
}; template <int len = 1> // Arbitrary sized bitset 16 8086598962041600 41472 8 3 2 2 1 1 1 1 1 1 1
void MonsterHunter(int root=1){ void solve(int n){ // solution using bitset<len> 17 74801040398884800 64512 6 3 2 2 1 1 1 1 1 1 1 1
set<Info> T(A+1, A+N+1); T.erase(A[root]); if(len < n){ solve<std::min(len*2, MAXLEN)>(n); return; } } 18 897612484786617600 103680 8 4 2 2 1 1 1 1 1 1 1 1
while(!T.empty()){
auto v = *T.rbegin(); T.erase(prev(T.end())); 6.15 Fast I/O, Fast Div, Fast Mod < 10^k prime # of prime < 10^k prime
int now = v.idx, nxt = Find(Par[v.idx]); // @TODO namespace io { // thanks to cgiosy -------------------------------------------------------------
UF[now] = nxt; T.erase(A[nxt]); A[nxt] += A[now]; const signed IS=1<<20; char I[IS+1],*J=I; 1 7 4 10 9999999967
if(nxt != root) T.insert(A[nxt]); inline void daer(){if(J>=I+IS-64){ 2 97 25 11 99999999977
} // @TODO char*p=I;do*p++=*J++; 3 997 168 12 999999999989
} while(J!=I+IS);p[read(0,p,I+IS-p)]=0;J=I;}} 4 9973 1229 13 9999999999971
template<int N=10,typename T=int>inline T getu(){ 5 99991 9592 14 99999999999973
6.12 Hook Length Formula daer();T x=0;int k=0;do x=x*10+*J-'0'; 6 999983 78498 15 999999999999989
int HookLength(const vector<int> &young){ while(*++J>='0'&&++k<N);++J;return x;} 7 9999991 664579 16 9999999999999937
if(young.empty()) return 1; template<int N=10,typename T=int>inline T geti(){ 8 99999989 5761455 17 99999999999999997
vector<int> len(young[0]); daer();bool e=*J=='-';J+=e;return(e?-1:1)*getu<N,T>();} 9 999999937 50847534 18 999999999999999989
ll num = 1, div = 1, cnt = 0; struct f{f(){I[read(0,I,IS)]=0;}}flu; };
for(int i=(int)young.size()-1; i>=0; i--){ struct FastMod{ // typedef __uint128_t L; 6.18 DLAS(Diversified Late Acceptance Search)
for(int j=0; j<young[i]; j++){ ull b, m; FastMod(ull b) : b(b), m(ull((L(1) << 64) / b)) {} template<class T, class U>
num = num * ++cnt % MOD; ull reduce(ull a){ // can be proven that 0 <= r < 2*b T incMod(T x, U mod) { x += 1; return x == mod ? 0 : x; }
div = div * (++len[j] + young[i] - j - 1) % MOD; ull q = (ull)((L(m) * a) >> 64), r = a - q * b; template<class Domain, class CoDomain, size_t LEN = 5>
} return r >= b ? r - b : r; pair<Domain, CoDomain> dlas( function<CoDomain(Domain&)> f,
} } }; function<void(Domain&)> mutate,
return num * Pow(div, MOD-2) % MOD; inline pair<uint32_t, uint32_t> Div(uint64_t a, uint32_t b){ Domain const& initial, u64 maxIdleIters = -1ULL) {
} if(__builtin_constant_p(b)) return {a/b, a%b}; array<Domain, 3> S{initial, initial, initial};
uint32_t lo=a, hi=a>>32; CoDomain curF = f(S[0]), minF = curF;
6.13 Floating Point Add (Kahan) __asm__("div %2" : "+a,a" (lo), "+d,d" (hi) : "r,m" (b)); size_t curPos = 0, minPos = 0, k = 0;
template<typename T> T float_sum(vector<T> v){ return {lo, hi}; // BOJ 27505, q r < 2^32 array<CoDomain, LEN> fitness; fitness.fill(curF);
T sum=0, c=0, y, t; // sort abs(v[i]) increase? } // divide 10M times in ~400ms for(u64 idleIters=0; idleIters<maxIdleIters; idleIters++){
for(T i:v) y=i-c, t=sum+y, c=(t-sum)-y, sum=t; ull mulmod(ull a, ull b, ull M){ // ~2x faster than int128 CoDomain prvF = curF;
return sum; //worst O(eps*N), avg O(eps*sqrtN) ll ret = a * b - M * ull(1.L / M * a * b); size_t newPos = incMod(curPos, 3);
}//dnc: worst O(eps*logN), avg O(eps*sqrtlogN) return ret + M * (ret < 0) - M * (ret >= (ll)M); if (newPos == minPos) newPos = incMod(newPos, 3);
} // safe for 0 ≤ a,b < M < (1<<63) when long double is 80bit Domain &curS = S[curPos], &newS = S[newPos];
6.14 Random, PBDS, Bit Trick, Bitset
newS = curS; mutate(newS); CoDomain newF = f(newS);
mt19937 6.16 DP Optimization if(newF < minF) idleIters=0, minPos=newPos, minF=newF;
rd((unsigned)chrono::steady_clock::now().time_since_epoch().count());
// Quadrangle Inequality : C(a, c)+C(b, d) ≤ C(a, d)+C(b, c) if(newF == curF || newF < *max_element(all(fitness))){
uniform_int_distribution<int> rnd_int(l, r); // rnd_int(rd)
// Monotonicity : C(b, c) ≤ C(a, d) curPos = newPos; curF = newF;
uniform_real_distribution<double> rnd_real(0, 1);// rnd_real(rd)
// CHT, DnC Opt(Quadrangle), Knuth(Quadrangle and Monotonicity) } CoDomain& fit = fitness[k]; k = incMod(k, LEN);
// ext/pb_ds/assoc_container.hpp, tree_policy.hpp, rope
// Knuth: K[i][j-1] <= K[i][j] <= K[i+1][j] if(curF > fit || curF < fit && curF < prvF) fit = curF;
// namespace __gnu_pbds (find_by_order, order_of_key)
// 1. Calculate D[i][i], K[i][i] } return { S[minPos], minF };
// namespace __gnu_cxx (append(str), substr(l, r), at(idx))
// 2. Calculate D[i][j], K[i][j] (i < j) } // 점수 최소화하는 함수, f: 상태의 점수를 반환
template <typename T> using ordered_set = tree<T, null_type,
// Another: D[i][j] = min(D[i-1][k] + C[k+1][j]), C quadrangle //dlas<state_type, score_type>(f, mutate, initial, maxIdleIter)
less<T>, rb_tree_tag, tree_order_statistics_node_update>;
// i=1..k j=n..1 k=K[i-1,j]..K[i,j+1] update, vnoi/icpc22_mn_c //initial:초기 상태, mutate:상태를 참조로 받아서 임의로 수정(반환X)
bool next_combination(T &bit, int N){
//maxIdleIters:지역 최적해에서 알마나 오래 기다릴지
T x = bit & -bit, y = bit + x;
bit = (((bit & ~y) / x) >> 1) | y; 6.17 Highly Composite Numbers, Large Prime
return (bit < (1LL << N)); } < 10^k number divisors 2 3 5 71113171923293137 6.19 Stack Hack
long long next_perm(long long v){ ------------------------------------------------------------- int main2(){ return 0; }
long long t = v | (v-1); 1 6 4 1 1 int main(){
return (t + 1) | (((~t & -~t) - 1) >> (__builtin_ctz(v) + 1)); 2 60 12 2 1 1 size_t sz = 1<<29; // 512MB
} // __builtin_clz/ctz/popcount 3 840 32 3 1 1 1 void* newstack = malloc(sz);
for(submask=mask; submask; submask=(submask-1)&mask); 4 7560 64 3 3 1 1 void* sp_dest = newstack + sz - sizeof(void*);
for(supermask=mask; supermask<(1<<n); 5 83160 128 3 3 1 1 1 asm __volatile__("movq %0, %%rax\n\t"
supermask=(supermask+1)|mask); 6 720720 240 4 2 1 1 1 1 "movq %%rsp , (%%rax)\n\t"
int frq(int n, int i) { int j, r = 0; // # of digit i in [1, n] 7 8648640 448 6 3 1 1 1 1 "movq %0, %%rsp\n\t": : "r"(sp_dest): );
for (j = 1; j <= n; j *= 10) if (n / j / 10 >= !i) r += (n / 8 73513440 768 5 3 1 1 1 1 1 main2();
10 / j - !i) * j + (n / j % 10 > i ? j : n / j % 10 == i ? n % 9 735134400 1344 6 3 2 1 1 1 1 asm __volatile__("pop %rsp\n\t");
j + 1 : 0); 10 6983776800 2304 5 3 2 1 1 1 1 1 return 0; }
Soongsil University – PS akgwi Page 24 of 25
√ √
6.20 Python Decimal • 삼각 치환: a2 − x2 에서는 x = a sin t, a2 + x2 에서는 x = a tan t, 7.5 Counting
from fractions import Fraction 범위 조심
from decimal import Decimal, getcontext • 입체 도형 부피: x = a, x = b 사이에서 단면의 넓이 함수 A(x)가 연속이면 조건 없음 단사 함수 전사 함수
부피는 ab A(x)dx kn
R
getcontext().prec = 250 # set precision 함수 일치 k!/(k − n)! k! × S2 (n, k)
N, two, itwo = 200, Decimal(2), Decimal(0.5) • (원주각법): 연속함수 f (x)가 [a, b]에서 f (x) ≥ 0일 때, f (x)와 두 직선 공 구별 X C(n + k − 1, n) C(k, n) C(n − 1, n − k) 조심
# sin(x) = sum (-1)^n x^(2n+1) / (2n+1)! x = a, x = b, 그리고 x축으로 둘러싸인 영역을 y축으로 회전시킨 부피는 상자 구별 X B(n, k) [n ≤ k] S2 (n, k)
Rb
# cos(x) = sum (-1)^n x^(2n) / (2n)! a 2πxf (x)dx 모두 구별 X Pk (n + k) [n ≤ k] Pk (n)
def angle(cosT): • 곡선 길이: f 0(x)가 [a, b]에서 연속이면, x = a, x = b사이의 곡선 길이는
Rbp
#given cos(theta) in decimal return theta 1 + [f 0(x)]2 dx
a
for i in range(N): cosT=((cosT+1)/two)**itwo
• 회전 곡면 넓이: x축으로 회전시킨 부피는 ab 2πf (x) 1 + [f 0(x)]2 dx
R p
sinT = (1-cosT*cosT)**itwo H R R ∂M ∂L
return sinT*(2**N) • C (Ldx + M dy) = D ( ∂x − ∂y )dxdy
pi = angle(Decimal(-1)) • where C is positively oriented, piecewise smooth, simple, closed; D • 단사 함수: 상자에 공 최대 1개, 전사 함수: 상자에 공 최소 1개
is the region inside C; L and M have continuous partial derivatives • 중복 조합: C(n + k − 1, n) = H(n, k)
in D. • 공 구별 X, 전사 함수에서 n = k = 0일 때 정답 1인 거 조심
6.21 Java I/O R • 교란 순열: 모든 i에 대해 π(i) 6= i가 성립하는 길이 n 순열 개수
// java.util.*, java.math.*, java.io.* f (x) f 0(x) f (x)dx 초항(0-based): 1,P0, 1, 2, 9, 44, 265, 1854
public class Main{ // BufferedReader, BufferedWriter sin x cos x − cos x 일반항: D(n) = n k
k=0 (−1) n!/k!, D(n) ≈ n!/e
public static void main(String[] args) throws IOException { cos x − sin x sin x 점화식: D(0) = 1; D(1) = 0; D(n) = (n − 1)(D(n − 1) + D(n − 2))
br=new BufferedReader(new InputStreamReader(System.in)); tan x sec2 x = 1 + tan2 x − ln | cos x| 생성함수(EGF): e−x /(1 − x)
bw=new BufferedWriter(new OutputStreamWriter(System.out)); csc x − csc x cot x ln | tan(x/2)| • 카탈란 수
String[] ar = br.readLine().split(" "); sec x sec x tan x ln | tan(x/2 + π/4)| 초항(0-based): 1, 1, 2, 5, 14, 42, 132, 429, 1430, 4862, 16796, 58786
int a=Integer.parseInt(ar[0]), b=Integer.parseInt(ar[1]); cot x − csc 2x ln | sin x|
√ √ 일반항: Cat(n) = Comb(2n, n)/(n P + 1) = C(2n, n) − C(2n, n + 1)
bw.write(String.valueOf(a+b)+"\n");br.close();bw.close(); arcsin x 1/ √ 1 − x2 x arcsin x + √1 − x2 점화식: Cat(0) = 1; Cat(n) = n−1 k=0 Cat(k) × Cat(n − k − 1)
ArrayList<Integer> a = new ArrayList<>(); √
arccos x −1/ 1 − x 2 x arccos x − 1 − x2 생성함수(OGF): (1 − 1 − 4x)/(2x)
a.add(1234); a.get(0); a.remove(a.size()-1); a.clear(); ln(x2 +1) 여는 괄호 n개, 닫는 괄호 k(≤ n)개 = C(n + k, k) × (n − k + 1)/(n + 1)
}} arctan x 1/(1 + x2 ) x arctan x −
√ 2 • 제 1종 스털링 수: k개의 사이클로 구성된 길이 n 순열 개수
csc−1 x −1/x 2
√ x −1 x csc−1 (x) + cosh−1 |x|
점화식: S1 (n, 0) = [n = 0]; S1 (n, k) = S1 (n − 1, k − 1) + S1 (n − 1, k) ×
7 Notes sec−1 x 2
1/x x − 1 x sec−1 (x) − cosh−1 |x| k
ln(x2 +1) P 생성함수(EGF): (− ln(1 − x)) /k!
(n − 1),
7.1 Triangles/Trigonometry cot−1 x −1/(1 + x2 ) x cot−1 x + 2
성질: n k=0 S1 (n, k) = n!
7.3 Zeta/Mobius Transform • 제 2종 스털링 수: n개의 원소를 P k개의 공집합이 아닌 부분집합으로 분할
• k차원 구 부피: V2k = π k /k!, v2k+1 = 2k+1 π k /(2k + 1)!!
일반항: S2 (n, k) = 1/k! × ki=1 (−1)k−i × C(k, i) × in , 단 S(0, 0) = 1
• sin(α ± β) = sin α cos β ± cos α sin β, sin 2θ = 2 sin θ cos θ • Subset Zeta/Mobius Transform(n=sz=2k ,i*=2)
점화식: S2 (n, 0) = [n = 0]; S2 (n, k) = S2 (n−1, k−1)+S2 (n−1, k)×k
• cos(α ± β) = cos α cos β ∓ sin α sin β, cos 2θ = 1 − 2 sin2 θ - for i=1..n-1 for j=0..n-1 if(i and j) v[j] ±= v[i xor j]
생성함수(EGF): (exp(x) − 1)k /k!
• tan α±tan β
tan(α ± β) = 1∓tan 2 tan θ
, tan 2θ = 1−tan • Superset Zeta/Mobius Transform(n=sz=2k ,i*=2)
α tan β 2θ 성질: A, B를 n−k, b(k−1)/2c의 켜진 비트 위치라고 하면, S2 (n, k) mod
- for i=1..n-1 for j=0..n-1 if(i and j) v[i xor j] ±= v[j]
• 반각 공식: sin2 θ/2 = 1−cos θ
, cos2 θ/2 = 1+cos θ
, tan2 θ/2 = 1−cos θ 2 = [A ∩ B = ∅]
2 2 1+cos θ • Divisor Zeta/Mobuis Transform(n=sz-1)
• 변 길이 a, p
b, c; p = (a + b + c)/2 성질: S2 (2n, n)은 n이 fibbinary number(연속한 1 없음) 일 때만 홀수
- for p:Prime for i=1..n/p v[i*p] += v[i]
• 넓이 A = p(p − a)(p − b)(p − c) • 벨 수 B(n): n개의 원소를 분할하는 방법의 수
- for p:Prime for i=n/p..1 v[i*p] -= v[i]
• 외접원 반지름 R = abc/4A, 내접원 반지름 r = A/p 초항(0-based): 1,P1, 2, 5, 15, 52, 203, 877, 4140, 21147, 115975
• Multiple Zeta/Mobius Transform(n=sz-1)
√ q
a 2 일반항: B(n) = n k=0 S2 (n, k), 몇 개의 상자를 버릴지 다 돌아보기
• 중선 길이 ma = 0.5 2b2 + 2c2 − a2 , 각 이등분선 sa = bc(1 − b+c ) - for p:Prime for i=n/p..1 v[i] += v[i*p]
점화식: B(0) = 1; B(n) = n−1
P
sinA - for p:Prime for i=1..n/p v[i] -= v[i*p] k=0 C(n − 1, k) × B(k)
• 사인 법칙 = 1/2R, 코사인 법칙 a2 = b2 + c2 − 2bc cos A 생성함수(EGF): exp(exp(x) − 1)
a • AND Convolution: SupZeta(A), SupZeta(B), SupMobius(mul)
• a+b
탄젠트 법칙 a−b
tan(A+B)/2
= tan(A−B)/2 • OR Conv.: Subset, GCD Conv.: Multiple, LCM Conv.: Divisor • 벨 수 B(n, k): n개를 집합 k개로 분할하는 방법의 수(공집합 허용)
n Pk−i (−1)j
• AND/OR: 22 0 0.3s, Subset: 220 2.5s, GCD/LCM: 1e6 0.3s 일반항: B(n, k) = ki=0 S2 (n, i) = ki=0 ii!
P P
• 중심 좌표 ( αxaα+β+γ
+βxb +γxc αya +βyb +γyc
, α+β+γ
) j=0 j!
• 분할 수 P (n): n을 자연수 몇 개의 합으로 나타내는 방법의 수, 순서 X
7.4 Generating Function 초항(0-based): 1, 1, 2, 3, 5, 7, 11, 15, 22, 30, 42 √
이름 α β γ • 등차수열: (pn + q)xn= p/(1 − x) + q/(1 − x)2 점화식: P (0) = 1, 팀노트 어딘가에 있는 코드로 O(n n) 가능
외심 a2 A b2 B c2 C A = b2 + c2 − a2 • 등비수열: (rx)n = (1 − rx)− 1 • 분할 수 P (n, k): n을 k개의 자연수의 합으로 나타내는 방법의 수, 순서 X
내심 a b c B = a2 + c2 − b2 • 조합: C(m, n)xn = (1 + x)m 점화식: PP (0, 0) = 1; P (n, 0) = 0; P (n, k) = P (n−1, k−1)+P (n−k, k)
무게중심 1 1 1 C = a2 + b2 − c2 • 중복조합: + n − 1, n)xn = (1 − x)−m
수심 BC CA AB PC(m
n
성질: n k=0 P (n, k) = P (n)
• f (n) = k=0 k! × S2 (n, k)의 EGF: 1/(2P− ex ) • 피보나치 수 F0 = 0, F1 = 1, Fn = Fn−1 + Fn−2
방심(A) −a b c • Ordinary Generating Function A(x) = i≥0 ai xi F =(φn −
일반항: n √ n √
7.2 Calculus, Newton’s Method A(rx) ⇒ rn an , xA(x)0 ⇒ nan φ )/ 5,φ={(1+ 5)/2} n
• 음함수 미분: f (x, y) = 0의 양변을 x에 대해 미분한 뒤 dy/dx에 대해 정리 A(x) + B(x) ⇒ an + bn , A(x)B(x) ⇒ n
P 생성함수(OGF): F (x) = x/(1 − x − x2 )
i=0 ai bn−i
• d
(예시) dx d
(x3 ) + dx d
(y 3 ) − 3 dx dy
(xy) = 3x2 + 3y 2 dx dy
− 3(y + x dx )=0 A(x)k ⇒ i1 +i2 +···+ik =n ai1 ai2 . . . aik
P 성질: F1 + · · · + Fn = Fn+2 − 1, F11 + · · · + Fn2 = Fn Fn+1
성질: gcd(Fn , Fm ) = Fgcd(n,m) , Fn2 − Fn+1 Fn−1 = (−1)n−1
• 역함수 미분: (f −1 )0(x) = 1/f 0(f −1 (x)) A(x)
⇒ n
P
1−x i=0 ai • 벨 수 B(n, k) 식 전개
• 뉴턴 랩슨: xn+1 = xn − f (xRn )/f 0 (xn ) R • Exponential Generating Function A(x) = i≥0 ai!i xi
P
B(n, k) = kj=0 S(n, j) = kj=0 1/j! ji=0 (−1)j−i jCi × in
P P P
• 치환 적분: x = g(t)로 두면, f (x)dx = f (g(t))g0(t)dt
A(x) + B(x) ⇒ an + bn , A(x)B(x) ⇒ i=0 ni ai bn−i
Pn
R f 0(x) Pk Pj (−1)j−i n Pk Pk (−1)j−i n
• (예시) dx에서 t = f (x)라고 두면 f 0(x) = dt/dx = j=0 i=0 i!(j−i)! i = i=0 j=i i!(j−i)! i
f (x) A(k) (x) ⇒Pan+k , xA(x) ⇒ nan
R f 0(x) n (−1)j n n Pk−i (−1)j
A(x)k ⇒ i +i +···+i =n = ki=0 k−i i = ki=0 i
R P P P
따라서 dx = 1/t dt = ln |t| + C = ln |f (x)| + C ai1 ai2 . . . aik j=0 j=0
Soongsil University – PS akgwi Page 25 of 25
Pn c
7.6 Faulhaber’s Formula ( k=1 k ) • 페르마 포인트 : 삼각형의 세 꼭짓점으로부터 거리의 합이 최소가 되는 점 7.11 About Graph Matching(Graph with |V | ≤ 500)
• Bn : 베르누이 수 2π/3 보다 큰 각이 있으면 그 점이 페르마 포인트, 그렇지 않으면 각 변마 • Game on a Graph : s에 토큰이 있음. 플레이어는 각자의 턴마다 토큰
1
• 생성함수(egf): exx−1 = (ex −1)/x = ∞ Bn n
P
n=0 n! x 다 정삼각형 그린 다음, 정삼각형의 끝점에서 반대쪽 삼각형의 꼭짓점으로 을 인접한 정점으로 옮기고 못 옮기면 짐.
Pn k!(−1)k 연결한 선분의 교점 s를 포함하지 않는 최대 매칭이 존재함 ↔ 후공이 이김
• 일반항: Bn = k=0 k+1 S2 (n, k) q √ • Chinese Postman Problem : 모든 간선을 방문하는 최소 가중치
1 Pn−1 n+1
2π/3 보다 큰 각이 없으면 거리의 합은 (a2 + b2 + c2 + 4 3S)/2, S
• 점화식: B0 = 1; Bn = − n+1 r=0 r
Br Walk를 구하는 문제.
는 넓이
Pn c =
Pc (−1)k c+1 c+1−k Floyd를 돌린 다음, 홀수 정점들을 모아서 최소 가중치 매칭 (홀수 정점은
• k B k n • 오일러 정리: 서로소인 두 정수 a, n에 대해 aφ(n) ≡ 1 (mod n)
k=1 k=0 c+1 k 짝수 개 존재)
모든 정수에 대해 an ≡ an−φ(n) (mod n)
7.7 About Graph Degree Sequence • Unweighted Edge Cover : 모든 정점을 덮는 가장 작은(minimum
m ≥ log2 n이면 am ≡ am%φ(n)+φ(n) (mod n)
• 단순 무향 그래프(Erdos-Gallai): 차수열 d1 ≥ · · · ≥ dn 의 합이 짝수 and cardinality/weight) 간선 집합을 구하는 문제
• g 0 + g 1 + g 2 + · · · g p−2 ≡ −1 (mod p) iff g = 1, otherwise 0.
|V | − |M |, 길이 3짜리 경로 없음, star graph 여러 개로 구성
모든 1 ≤ k ≤ n에 대해 ki=1 di ≤ k(k − 1) + n
P P
i=k+1 min(di , k) • if n ≡ 0 (mod 2), then 1n + 2n + · · · + (n − 1)n ≡ 0 (mod n)
• Weighted Edge Cover : sumv∈V (w(v)) − sum(u,v)∈M (w(u) +
• 단순 이분 그래프(Gale-Ryser): 차수열 a1 ≥ · · · ≥ an , bi 에서 sum(a) = • Eulerian numbers
Pk Pn w(v) − d(u, v)), w(x)는 x와 인접한 간선의 최소 가중치
sum(b) and 모든 1 ≤ k ≤ n에 대해 i=1 ai ≤ i=1 min(bi , k) Number of permutations π ∈ Sn in which exactly k elements are
• NEERC’18 B : 각 기계마다 2명의 노동자가 다뤄야 하는 문제.
• 단순 방향 그래프(Fulkerson-Chen-Anstee): a1 ≥ · · · ≥ an 를 만족 greater than the previous element. k j:s s.t. π(j) > π(j + 1), k + 1
기계마다 두 개의 정점을 만들고 간선으로 연결하면 정답은 |M | − |기계|
하는 진입/진출 차수열 (a1 , b1 ), · · · , (an , bn )에서 sum(a) = sum(b) j:s s.t. π(j) ≥ j, k j:s s.t. π(j) > j.
Pk Pk 임. 정답에 1/2씩 기여한다는 점을 생각해보면 좋음.
E(n, k) = (n − k)E(n − 1, k − 1) + (k + 1)E(n − 1, k)
Pn 모든 1 ≤ k ≤ n에 대해
and i=1 ai ≤ i=1 min(bi , k − 1) + • Min Disjoint Cycle Cover : 정점이 중복되지 않으면서 모든 정점을
i=k+1 min(bi , k)
E(n, 0) = E(n, n − 1) = 1
덮는 길이 3 이상의 사이클 집합을 찾는 문제.
E(n, k) = kj=0 (−1)j n+1 (k + 1 − j)n
P
7.8 Burnside, Grundy, Pick, Hall, Simpson, Area of j 모든 정점은 2개의 서로 다른 간선, 일부 간선은 양쪽 끝점과 매칭되어야
• Pythagorean triple: a + b = c2 이고 서로소인 (a, b, c) 생성
2 2
하므로 플로우를 생각할 수 있지만 용량 2짜리 간선에 유량을 1만큼 흘릴
Quadrangle, Fermat Point, Euler, Pythagorean s2 −t2 s2 +t2 수 있으므로 플로우는 불가능.
• Burnside’s Lemma (a, b, c) = (st, 2
, 2 ), gcd(s, t) = 1, s > t
각 정점과 간선을 2개씩((v, v 0 ), (ei,u , ei,v ))로 복사하자. 모든 간선
- 수식
e = (u, v)에 대해 eu 와 ev 를 잇는 가중치 w짜리 간선을 만들고(like
G=(X,A): 집합X와 액션A로 정의되는 군G에 대해, |A||X/A| = 7.9 About Graph Minimum Cut NEERC18), (u, ei,u ), (u0 , ei,u ), (v, ei,v ), (v 0 , ei,v )를 연결하는 가중치 0
sum(|Fixed points of a|, for all a in A) • N 개의 boolean 변수 v1 , · · · , vn 을 정해서 비용을 최소화하는 문제 짜리 간선을 만들자. Perfect 매칭이 존재함 ⇔ Disjoint Cycle Cover 존재.
X/A 는 Action으로 서로 변형가능한 X의 원소들을 동치로 묶었을때 동 =true인 점은 T , false인 점은 F 와 연결되게 분할하는 민컷 문제 최대 가중치 매칭 찾은 뒤 모든 간선 가중치 합에서 매칭 빼면 됨.
치류(파티션) 집합
• Two Matching : 각 정점이 최대 2개의 간선과 인접할 수 있는 최대
- 풀어쓰기 1. vi 가 T일 때 비용 발생: i에서 F로 가는 비용 간선 가중치 매칭 문제.
orbit: 그룹에 대해 두 원소 a,b와 액션f에 대해 f(a)=b인거에 간선연결한 2. vi 가 F일 때 비용 발생: i에서 T로 가는 비용 간선 각 컴포넌트는 정점 하나/경로/사이클이 되어야 함. 모든 서로 다른 정점
컴포넌트(연결집합) 3. vi 가 T이고 vj 가 F일 때 비용 발생: i에서 j로 가는 비용 간선 쌍에 대해 가중치 0짜리 간선 만들고, 가중치 0짜리 (v, v 0 ) 간선 만들면
orbit개수 = sum(각 액션 g에 대해 f(x)=x인 x(고정점)개수)/액션개수 4. vi 6= vj 일 때 비용 발생: i에서 j로, j에서 i로 가는 비용 간선 Disjoing Cycle Cover 문제가 됨. 정점 하나만 있는 컴포넌트는 self-loop,
- 자유도 치트시트 5. vi 가 T면 vj 도 T여야 함: i에서 j로 가는 무한 간선 경로 형태의 컴포넌트는 양쪽 끝점을 연결한다고 생각하면 편함.
회전 n개: 회전i의 고정점 자유도=gcd(n,i) 6. vi 가 F면 vj 도 F여야 함: j에서 i로 가는 무한 간선
임의뒤집기 n=홀수: n개 원소중심축(자유도 (n+1)/2)
• 5/6번 + vi 와 vj 가 달라야 한다는 조건이 있으면 MAX-2SAT 7.12 Checklist
임의뒤집기 n=짝수: n/2개 원소중심축(자유도 n/2+1) + n/2개 원소안
• Maximum Density Subgraph (NEERC’06H, BOJ 3611 팀의 난이도) • (예비소집) bits/stdc++.h, int128, long double 80bit, avx2 확인
지나는축(자유도 n/2)
• (예비소집) 스택 메모리(지역 변수, 재귀, 람다 재귀), 제출 파일 크기 확인
• 알고리즘 게임
– density ≥ x인 subgraph가 있는지 이분 탐색 • (예비소집) MLE(힙,스택), stderr 출력 RTE?, 줄 앞뒤 공백 채점 결과
- Nim Game의 해법(마지막에 가져가는 사람이 승) : XOR = 0이면 후공
– 정점 N 개, 간선 M 개, 차수 Di 개 • 비슷한 문제를 풀어본 적이 있는가?
승, 0 아니면 선공 승
– 그래프의 간선마다 용량 1인 양방향 간선 추가 • 단순한 방법에서 시작할 수 있을까? (Brute Force)
- Subtraction Game : 한 번에 k 개까지의 돌만 가져갈 수 있는 경우, 각
– 소스에서 정점으로 용량 M , 정점에서 싱크로 용량 M − Di + 2x • 내가 문제를 푸는 과정을 수식화할 수 있을까? (예제를 직접 해결하면서)
더미의 돌의 개수를 k + 1로 나눈 나머지를 XOR 합하여 판단한다.
– min cut에서 S와 붙어 있는 애들이 1개 이상이면 x 이상이고, 그게 • 문제를 단순화할 수 없을까? / 그림으로 그려볼 수 있을까?
- Index-k Nim : 한 번에 최대 k개의 더미를 골라 각각의 더미에서 아무
subgraph의 정점들 • 수식으로 표현할 수 있을까? / 문제를 분해할 수 있을까?
렇게나 돌을 제거할 수 있을 때, 각 binary digit에 대하여 합을 k + 1로
– while(r-l ≥ 1.0/(n*n)) 으로 해야 함. 너무 많이 돌리면 실수 오차 • 뒤에서부터 생각해서 풀 수 있을까? / 순서를 강제할 수 있을까?
나눈 나머지를 계산한다. 만약 이 나머지가 모든 digit에 대하여 0이라면
• 특정 형태의 답만을 고려할 수 있을까? (정규화)
두번째, 하나라도 0이 아니라면 첫번째 플레이어가 승리.
• 구간을 통째로 가져간다 : 플로우 + 적당한 자료구조
- Misere Nim : 모든 돌 무더기가 1이면 N이 홀수일 때 후공 승, 그렇지 7.10 Matrix with Graph(Kirchhoff, Tutte, LGV) (i, i + 1, k, 0), (s, e, 1, w), (N, T, k, 0)
않은 경우 XOR 합 0이면 후공 승 • Kirchhoff’s Theorem : 그래프의 스패닝 트리 개수 • a = b : a만 이동, b만 이동, 두 개 동시에 이동, 반대로 이동
• Pick’s Theorem - m[i][j] := -(i-j 간선 개수) (i ≠ j) / (유향) -(i → j 간선) • 말도 안 되는 것 / 당연하다고 생각한 것 다시 생각해 보기
격자점으로 구성된 simple polygon이 주어짐. I 는 polygon 내부의 격자 - m[i][i] := 정점 i의 degree / (유향) 정점 i의 in degree • Directed MST / Dominator Tree
점 수, B 는 polygon 선분 위 격자점 수, A는 polygon의 넓이라고 할 때, - res = (m의 첫 번째 행과 첫 번째 열을 없앤 (n-1) by (n-1) matrix의 • 일정 비율 충족 or 2 3개로 모두 커버 : 랜덤
다음과 같은 식이 성립한다. A = I + B/2 − 1 행렬식) • 확률 : DP, 이분 탐색(NYPC 2019 Finals C)
• 홀의 결혼 정리 : 이분그래프(L-R)에서, 모든 L을 매칭하는 필요충분 조건 - (유향) m의 루트 번째 행과 열을 삭제한 행렬의 행렬식 • 최대/최소 : 이분 탐색, 그리디(Prefix 고정, Exchange Argument),
= L에서 임의의 부분집합 S를 골랐을 때, 반드시 (S의 크기) <= (S와 • Tutte Matrix : 그래프의 최대 매칭 DP(순서 고정)
연결되어있는 모든 R의 크기)이다. - m[i][j] := 간선 (i, j)가 없으면 0, 있으면 i < j?r : −r, r은 [0, P ) 구간의
h • 냅색: 파라미터 순서 변경, min plus convolution, FFT
• Simpson 공식 (적분) P : Simpson 공식, Sn (f ) = 3 [f (x0 ) + f (xn ) + 임의의 정수
P • signal(SIGSEGV,[](int){ Exit(0);}); converts segfaults into WA.
4 f (x2i+1 ) + 2 f (x2i )] - rank(m)/2가 높은 확률로 최대 매칭 (mod P)
M (b−a)
SIGABRT(assertion fail), SIGFPE(0div)
- M = max |f 4 (x)|이라고 하면 오차 범위는 최대 En ≤ 180 h4 • LGV Theorem: 간선에 가중치 있는 DAG에서 어떤 경로 P 의 간선 가중 • feenableexcept(29) kills problem on NaNs(1), 0div(4), inf(8), denor-
• 브라마굽타 : 원에 내접하는 p 사각형의 각 선분의 길이가 a, b, c, d일 때 치 곱을 w(P ), 모든 a → b 경로들의 w(P )의 합을 e(a, b)라고 하자. n mals(16)
사각형의 넓이 S = (s − a)(s − b)(s − c)(s − d), s = (a + b + c + d)/2 개의 시작점 ai 와 도착점 bj 가 주어졌을 때, 서로 정점이 겹치지 않는 n
• 브레치나이더 : 임의의 사각형의 각 변의 길이를 a, b, c, d라고 하 개의 경로로 시작점과 도착점을 일대일 대응시키는 모든 경우에서 w(P )
고,
p 마주보는 두 각의 합을 2로 나눈 값을 θ라 하면, S = 의 곱의 합은 det M (i, j) = e(ai , bj )와 같음. 따라서 모든 가중치를 1로
(s − a)(s − b)(s − c)(s − d) − abcd × cos2 θ 두면 서로소 경로 경우의 수를 구함