0.0.1

cg/polygon.cpp

@@ -0,0 +1,81 @@
+#define nxt(x) ((x + 1) >= n ? x + 1 - n : x + 1)
+struct polygon {    
+    vector<vec> s;
+    polygon () {s.clear();}
+    ld circ () const { // 多边形的周长
+        if (s.size() <= 1) return 0.0;
+        ld ans = 0.0;
+        for (int i = 0, n = s.size(); i < n; ++i) 
+            ans += s[i].dis(s[nxt(i)]);
+        return ans;
+    }
+    ld area () const { // 多边形的面积
+        if (s.size() <= 1) return 0.0;
+        ld ans = 0.0;
+        for (int i = 0, n = s.size(); i < n; ++i)
+            ans += s[i].cross(s[nxt(i)]);
+        return ans / 2;
+    }
+    void polar_angle_sort_origin() { // 关于原点极角排序
+        auto cmp = [&](const vec &a, const vec &b) {
+            if (a.up() != b.up()) return a.up() > b.up(); 
+            ld d = a.cross(b); // 叉积判断方向
+            return d == 0 ? a.norm() < b.norm() : d > 0;
+        };
+        sort(s.begin(), s.end(), cmp);
+    }
+    void polar_angle_sort() { // 极角排序
+        auto cmp = [&](const vec &a, const vec &b) {
+            double d = (a - s[0]).cross(b - s[0]);
+            return d == 0 ? a.dis2(s[0]) < b.dis2(s[0]) : d > 0;
+        };
+        for (int i = 1, n = s.size(); i < n; ++i)
+            if (s[i] < s[0]) swap(s[0], s[i]);
+        sort(s.begin() + 1, s.end(), cmp); // 注意不要把 s[0] 排进去
+    }
+    void convex_hull () { // 求凸包
+        if (s.size() <= 1) return;
+        polar_angle_sort(); // 先进行极角排序
+        vector<vec> t; t.clear(); swap(s, t); 
+        s.push_back(t[0]); s.push_back(t[1]); 
+        for (int top = 1, i = 2, n = t.size(); i < n; ++i) {
+            while ((s[top] - s[top - 1]).cross(t[i] - s[top - 1]) <= 0) {
+                s.pop_back(); if ((--top) < 1) break;
+            }
+            s.push_back(t[i]); ++top;
+        }
+    }
+    void convex_hull_maxcnt() { // 上下凸壳合并版写法
+        auto cmp = [&](const vec &a, const vec &b) {
+            return a.y == b.y ? a.x < b.x : a.y < b.y;
+        };
+        sort(s.begin(), s.end(), cmp);
+        vector<vec> t; t.clear(); swap(s, t);
+        s.push_back(t[0]); s.push_back(t[1]);
+        int top = 1;
+        for (int i = 2, n = t.size(); i < n; ++i) { // 保留最多的点
+            while ((s[top] - s[top - 1]).cross(t[i] - s[top - 1]) < 0) {
+                s.pop_back(); if ((--top) < 1) break;
+            }
+            s.push_back(t[i]); ++top;
+        }
+        int lim = top;
+        for (int i = t.size() - 2; i >= 0; --i) {
+            while ((s[top] - s[top - 1]).cross(t[i] - s[top - 1]) < 0) {
+                s.pop_back(); if ((--top) < lim) break;
+            }
+            s.push_back(t[i]); ++top;
+        }
+        s.pop_back();
+    }
+    ld convex_diam () const { // 凸包的直径
+        if (s.size() <= 1) return 0.0;
+        ld ans = 0;
+        for (int i = 0, ptr = 1, n = s.size(); i < n; ++i) {
+            while (s[i].dis(s[ptr]) <= s[i].dis(s[nxt(ptr)])) 
+                ptr = nxt(ptr);
+            ans = max(ans, s[i].dis(s[ptr]));
+        }
+        return ans;
+    }
+};

cg/seg.cpp

@@ -0,0 +1,32 @@
+struct seg {
+    vec p1, p2;
+    vec proj(const vec &b) {return p1 + (p2 - p1).proj(b - p1);}
+    vec refl(const vec &b) {return p1 + (p2 - p1).refl(b - p1);}
+    bool cross(const seg &b) const{ // 判断两线段是否相交
+        ld w1 = (p2 - p1).cross(b.p1 - p1);
+        ld w2 = (p2 - p1).cross(b.p2 - p1);
+        ld w3 = (b.p2 - b.p1).cross(p1 - b.p1);
+        ld w4 = (b.p2 - b.p1).cross(p2 - b.p1);
+        if (w1 * w2 > 0 || w3 * w4 > 0) return 0;
+        if (w1 * w2 < 0 || w3 * w4 < 0) return 1;
+        ld v1 = (b.p1 - p1).dot(b.p1 - p2);
+        ld v2 = (b.p2 - p1).dot(b.p2 - p2);
+        ld v3 = (p1 - b.p1).dot(p1 - b.p2);
+        ld v4 = (p2 - b.p1).dot(p2 - b.p2);
+        return (v1 <= 0 || v2 <= 0 || v3 <= 0 || v4 <= 0);
+    }
+    vec cross_point(const seg &b) const { // 相交且不平行的情况下求交点
+        vec v = p2 - p1, u1 = b.p1 - p1, u2 = b.p2 - p1;
+        ld s1 = fabs(v.cross(u1)) / 2, s2 = fabs(v.cross(u2)) / 2;
+        return p1 + u1 * (s2 / (s1 + s2)) + u2 * (s1 / (s1 + s2));
+    }
+    ld dis(const vec &b) const { // 点到线段距离
+        if ((p2 - p1).dot(b - p1) < 0 || (p1 - p2).dot(b - p2) < 0)
+            return min(p1.dis(b), p2.dis(b)); // 垂足不在线段上
+        return fabs((p2 - p1).cross(b - p1)) / (p2 - p1).len(); 
+    }
+    ld dis(const seg &b) const { // 线段到线段距离
+        if (cross(b)) return 0.0;
+        return min(min(dis(b.p1), dis(b.p2)), min(b.dis(p1), b.dis(p2)));
+    }
+};

cg/vec.cpp

@@ -0,0 +1,31 @@
+const ld a90 = atan(1) * 2;
+const ld a180 = a90 * 2, a270 = a90 * 3, a360 = a90 * 4;
+const ld dlt[2][2] = {0, a180, a360, a180};
+struct vec {
+    ld x, y;
+    vec (ld _x = 0, ld _y = 0){x = _x; y = _y;}
+    vec operator + (const vec &b) const {return (vec){x + b.x, y + b.y};} 
+    vec operator - (const vec &b) const {return (vec){x - b.x, y - b.y};} 
+    vec operator * (const ld &k) const {return (vec){x * k, y * k};}
+    ld dot (const vec &b) const {return x * b.x + y * b.y;} 
+    ld cross (const vec &b) const {return x * b.y - y * b.x;}
+    ld norm () const {return x * x + y * y;}
+    ld len () const {return sqrt(norm());}
+    ld dis (const vec &b) const {return ((*this) - b).len();}
+    ld dis2 (const vec &b) const {return ((*this) - b).norm();}    
+    ld angle() const { // 与 x 正半轴夹角，逆时针旋转
+        if (fabs(x) < 1e-15) return a90 * (1 + 2 * (y < 0));
+        else return atan(y / x) + dlt[y < 0][x < 0];
+    }
+    vec turn (const ld &ang) { // 逆时针旋转 ang 角度
+        ld cosa = cos(ang), sina = sin(ang);
+        return (vec){x * cosa - y * sina, x * sina + y * cosa};
+    }
+    vec proj (const vec &b) {return (*this) * (dot(b) / this -> norm());}
+    vec refl (const vec &b) {return proj(b) * 2 - b;} //对称
+    //找最下的点里，最靠左的，把其他点控制在一二象限和x正半轴
+    bool up() const {return (y > 0 || (y == 0 && x > 0));}
+    bool operator < (const vec &b) const{
+        return y == b.y ? x < b.x : y < b.y;
+    }
+};

graph/scc.cpp

@@ -0,0 +1,40 @@
+bool vis[N];
+int n, m, scc, bl[N];
+vector<int> e[N], re[N], nodes[N], order;
+void add(int u, int v) {e[u].pb(v); re[v].pb(u);}
+void dfs1(int u) {
+    vis[u] = 1;
+    for (auto v : e[u]) if (!vis[v]) dfs1(v);
+    order.pb(u);
+}
+void dfs2(int u) {
+    nodes[bl[u] = scc].pb(u);
+    for (auto v : re[u]) if (!bl[v]) dfs2(v);
+}
+void kosaraju(int n) {
+    for (int u = 1; u <= n; ++u) if (!vis[u]) dfs1(u);
+    reverse(order.begin(), order.end());
+    for (auto u : order) if (!bl[u]) {++scc; dfs2(u);}
+}
+// ===== shrink nodes =====
+int indeg[N], outdeg[N];
+vector<int> source, sink;
+unordered_map<ll, bool> valid;
+inline bool has_edge(int u, int v) {
+    ll w = 1000000000ll * u + v;
+    if (valid[w]) return 1; 
+    valid[w] = 1; return 0;
+}
+inline void shrink(int n) {
+    valid.clear();
+    for (int u = 1; u <= n; ++u) e[u].clear();
+    for (int v = 1; v <= n; ++v)
+        for (auto u : re[v]) 
+            if (bl[u] != bl[v] && !has_edge(bl[u], bl[v])) {
+                ++indeg[bl[v]]; ++outdeg[bl[u]]; e[bl[u]].pb(bl[v]);
+            }
+    for (int u = 1; u <= scc; ++u) if (!indeg[u]) source.pb(u); 
+    for (int u = 1; u <= scc; ++u) if (!outdeg[u]) sink.pb(u);
+    for (int u = 1; u <= n; ++u) re[u].clear();
+    for (int u = 1; u <= scc; ++u) for (auto v : e[u]) re[v].pb(u);
+}

math/comb.cpp

@@ -0,0 +1,22 @@
+namespace Comb {
+    int fac[N], ifac[N];
+
+    inline int fpow(int x, int t) {
+        int res = 1;
+        for (; t; t >>= 1, x = 1ll * x * x % mod)
+            if (t & 1) res = 1ll * res * x % mod;
+        return res;
+    }
+
+    inline void init() {
+        fac[0] = ifac[0] = 1;
+        for (int i = 1; i < N; ++i) fac[i] = 1ll * fac[i - 1] * i % mod;
+        ifac[N - 1] = fpow(fac[N - 1], mod - 2);
+        for (int i = N - 2; i; --i) ifac[i] = 1ll * ifac[i + 1] * (i + 1) % mod;
+    }
+
+    inline int C(int n, int m) {
+        if (n < m) return 0;
+        return 1ll * fac[n] * ifac[m] % mod * ifac[n - m] % mod;
+    }
+}

math/linear_sieve.cpp

@@ -0,0 +1,22 @@
+int mnd[N], mncnt[N], mnpw[N], prm[N], phi[N] = {0, 1};
+for (int i = 2; i < N; ++i) {
+    if (!mind[i]) {
+        prm[++prm[0]] = mnd[i] = mnpw[i] = i; mncnt[i] = 1;  
+        phi[i] = i - 1; // cntd[i] = 2; // sumd[i] = 1 + i;  
+    }
+    for (int j = 1, p = prm[1], n; j <= prm[0]; p = prm[++j]) {
+        if ((n = i * p) > N) break;
+        if (p == mind[i]) {
+            mind[n] = p;
+            pw[n] = pw[i] * p; 
+            cnt[n] = cnt[i] + 1;
+            if (i == pw[i]) phi[n] = phi[i] * p; 
+                            // cntd[n] = cnt[n] + 1; 
+                            // sumd[n] = (pw[n] * p - 1) / (p - 1);
+            else phi[n] = phi[i / pw[i]] * phi[pw[i] * p];
+            break;
+        }
+        mnd[n] = pw[n] = p; cnt[n] = 1; 
+        phi[n] = phi[i] * phi[p];
+    }
+}