CF#538 E. Arithmetic Progression

問題 https://codeforces.com/contest/1114/problem/E
解説 https://codeforces.com/blog/entry/65136


初項 a[0] は わかる。0 から n-1 の乱数 r[0],...,r[29] を つくって d=gcd(a[r[0]] - a[0], a[r[1]] - a[0], ... , a[r[29]] - a[0]) をもとめると 高確率で それがこたえである。確率をしらべよう。こたえがまちがっているとき r[0] から r[29] はすべて ある倍数になっている。r[0] から r[29] がすべて d の倍数であるという事象を A[d] とすると まちがえる確率は

\begin{align}
P(A_2 \cup A_3 \cup A_4 \cup A_5 \cup \cdots ) &\le P(A_2)+P(A_3)+P(A_4)+P(A_5)+\cdots \\
&\le \frac{1}{2^{29}} + \frac{1}{3^{29}} + \frac{1}{4^{29}} + \frac{1}{5^{29}} + \cdots \\
&= \zeta(29)-1 \\
&\le 1.8627 \times 10^{-9}
\end{align}

で おさえられる。解析のしやすさを かんがえて 乱数で初項は あらわれないと仮定して計算している。30 ではなく 29 にしているのは この仮定をとりのぞくため。積分でおさえてもいいかも というか 積分の方が てで やるには むいている。

\[ \bigcirc \le \frac{1}{2^{29}} + \int_{2}^{\infty} \frac{dx}{x^{29}} \le \frac{1}{2^{28}} \]



Codeforces の乱数は わなが おおいです。


さて

1 から 1000 のなかから 重複なしで 10 個えらぶのと 重複ありで 10 個えらぶのと すべてが偶数になる確率がたかいのは どちらでしょう。

CF#537 | Destroy the Colony

https://codeforces.com/contest/1111/problem/D

ならびは あとであたえればいいから ただのナップサック。ならびを おぼえたままでも DP できるけど テイスウバイでしんだから えだがりして とおした。でも ナップサックのホウが ラクだね。ハイレツをつかいまわす DP でも もどせて そっちのホウが かきやすそうだし そうすればよかった。

おもしろいとおもった かんがえかたとして ナップサックのかぞえあげは

\begin{align}
(1+x^\bigcirc)(1+x^\bigcirc) \cdots (1+x^\bigcirc)
\end{align}

というボカンスウを つかうやりかたがあるけど もどすというのが $(1+x^\bigcirc)^{-1}$ をかけることに むすびつく。yukicoder の wiki にもかいてあるけど。

#include <stdio.h>
#include <string.h>

#define MOD 1000000007

struct mint {
    int n;
    mint(int n_ = 0) : n(n_) {}
};

mint operator+(mint a, mint b) { a.n += b.n; if (a.n >= MOD) a.n -= MOD; return a; }
mint operator-(mint a, mint b) { a.n -= b.n; if (a.n < 0) a.n += MOD; return a; }
mint operator*(mint a, mint b) { return (long long)a.n * b.n % MOD; }
mint &operator+=(mint &a, mint b) { return a = a + b; }
mint &operator-=(mint &a, mint b) { return a = a - b; }
mint &operator*=(mint &a, mint b) { return a = a * b; }

mint F[100001] = {1, 1};
mint IF[100001] = {1, 1};
mint I[100001] = {0, 1};

mint dp[53][50001];
int N;
char S[100001];
int A[52];
int A_[52];

int min(int x, int y) {
    return x < y ? x : y;
}

int f(int c) {
    if ('a' <= c && c <= 'z') {
        return c - 'a';
    } else {
        return c - 'A' + 26;
    }
}

void susume(int i) {
    if (A[i] == 0) {
        for (int j = 0; j <= N / 2; j++) {
            dp[i + 1][j] = dp[i][j];
        }
    } else {
        for (int j = 0; j <= N / 2; j++) {
            if (j >= A[i]) {
                dp[i + 1][j] = dp[i][j] + dp[i][j - A[i]];
            } else {
                dp[i + 1][j] = dp[i][j];
            }
        }
    }
}

void modore(int i) {
    for (int j = 0; j <= N / 2; j++) {
        if (j >= A[i]) {
            dp[i][j] = dp[i + 1][j] - dp[i][j - A[i]];
        } else {
            dp[i][j] = dp[i + 1][j];
        }
    }
}

int main() {
    for (int i = 2; i <= 100000; i++) {
        I[i] = I[MOD % i] * (MOD - MOD / i);
        F[i] = i * F[i - 1];
        IF[i] = I[i] * IF[i - 1];
    }
    scanf("%s", S);
    N = strlen(S);
    for (int i = 0; i < N; i++) {
        A[f(S[i])]++;
    }
    dp[0][0] = 1;
    for (int i = 0; i < 52; i++) {
        susume(i);
    }
    for (int i = 0; i < 52; i++) {
        A_[i] = A[i];
    }
    static mint ans[52][52];
    for (int x = 0; x < 52; x++) {
        if (A_[x] == 0) continue;
        ans[x][x] = dp[52][N / 2];
        A[51] = A_[x];
        modore(51);
        for (int y = x + 1; y < 52; y++) {
            if (A_[y] == 0) continue;
            A[50] = A_[y];
            modore(50);
            int s = A_[x] + A_[y];
            if (s > N / 2) {
                ans[x][y] = ans[y][x] = 0;
            } else {
                ans[x][y] = ans[y][x] = dp[50][N / 2] + dp[50][N / 2 - s];
            }
        }
    }
    mint v = F[N / 2] * F[N / 2];
    for (int i = 0; i < 52; i++) {
        v *= IF[A_[i]];
    }
    int q;
    scanf("%d", &q);
    while (q--) {
        int x, y;
        scanf("%d %d", &x, &y);
        x--; y--;
        x = f(S[x]);
        y = f(S[y]);
        printf("%d\n", (v * ans[x][y]).n);
    }
    return 0;
}

ケイシキテキベキキュウスウ そのばで できるセット。

$(1+x)$ と $(1+x)^{-1}$

#include <stdio.h>
#include <assert.h>

#define N 20
int a[N];

int main() {
    for (int i = 0; i < N; i++) {
        a[i] = i;
    }
    /* kakeru (x+1) */
    for (int i = N - 2; i >= 0; i--) {
        a[i + 1] += a[i];
    }
    /* waru (x+1) */
    for (int i = 0; i + 1 < N; i++) {
        a[i + 1] -= a[i];
    }
    for (int i = 0; i < N; i++) {
        assert(a[i] == i);
    }
    return 0;
}

$(1+x+x^2+\cdots)$ と $(1+x+x^2+\cdots)^{-1}$

#include <stdio.h>
#include <assert.h>

#define N 20
int a[N];

int main() {
    for (int i = 0; i < N; i++) {
        a[i] = i;
    }
    /* kakeru (1+x+x^2+x^3...) */
    for (int i = 0; i + 1 < N; i++) {
        a[i + 1] += a[i];
    }
    /* waru (1+x+x^2+x^3+...) */
    for (int i = N - 2; i >= 0; i--) {
        a[i + 1] -= a[i];
    }
    for (int i = 0; i < N; i++) {
        assert(a[i] == i);
    }
    return 0;
}

たとえば (1+x)^5 を かきだしたあと (1+x)^4 をかきだす。

#include <stdio.h>
#include <assert.h>

#define N 20
int a[N];

int main() {
    a[0] = 1;
    for (int j = 0; j < 5; j++) {
        for (int i = N - 2; i >= 0; i--) {
            a[i + 1] += a[i];
        }
    }
    for (int i = 0; i < N; i++) {
        printf("%d ", a[i]);
    }
    putchar('\n');
    for (int i = 0; i + 1 < N; i++) {
        a[i + 1] -= a[i];
    }
    for (int i = 0; i < N; i++) {
        printf("%d ", a[i]);
    }
    putchar('\n');
    return 0;
}
1 5 10 10 5 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 4 6 4 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

ついで ケイシキテキディリクレキュウスウ そのばで できるセット。

$\zeta(s)$ と $\mu(s)$

#include <stdio.h>
#include <assert.h>

#define N 20
int a[N];

int main() {
    for (int i = 1; i < N; i++) {
        a[i] = i;
    }
    /* kakeru zeta */
    for (int i = N - 1; i >= 1; i--) {
        for (int j = i * 2; j < N; j += i) {
            a[j] += a[i];
        }
    }
    /* waru zeta = kakeru moebius */
    for (int i = 1; i < N; i++) {
        for (int j = i * 2; j < N; j += i) {
            a[j] -= a[i];
        }
    }
    for (int i = 1; i < N; i++) {
        assert(a[i] == i);
    }
    return 0;
}

たとえば オイラーの phi は $\zeta(s-1)\mu(s)$ なので

#include <stdio.h>
#include <assert.h>

#define N 20
int a[N];

int main() {
    for (int i = 1; i < N; i++) {
        a[i] = i;
    }
    for (int i = 1; i < N; i++) {
        for (int j = i * 2; j < N; j += i) {
            a[j] -= a[i];
        }
    }
    for (int i = 1; i < N; i++) {
        printf("%d ", a[i]);
    }
    putchar('\n');
    return 0;
}
1 1 2 2 4 2 6 4 6 4 10 4 12 6 8 8 16 6 18

CF#537 | Tree

https://codeforces.com/contest/1111/problem/E

ねに ちかいノードから グループをきめていく。ノードたちを ねからの とおさで ならびかえて そのならびで DP する。ジョウタイは はじめの i コのノードのグループをきめた j コのグループがある のふたつ。センイは あたらしいグループをつくる いまあるグループにいれる のふたつ。ふたつめのセンイだけど おやノードのかずを k としたとき その k グループには いれられないから j-k とおりの えらびかたがある。

ひさびさに HL ブンカイかいたな。

#include <iostream>
#include <algorithm>
#include <vector>

using namespace std;

#define MOD 1000000007

struct mint {
    int n;
    mint(int n_ = 0) : n(n_) {}
};

mint operator+(mint a, mint b) { a.n += b.n; if (a.n >= MOD) a.n -= MOD; return a; }
mint operator-(mint a, mint b) { a.n -= b.n; if (a.n < 0) a.n += MOD; return a; }
mint operator*(mint a, mint b) { return (long long)a.n * b.n % MOD; }
mint &operator+=(mint &a, mint b) { return a = a + b; }
mint &operator-=(mint &a, mint b) { return a = a - b; }
mint &operator*=(mint &a, mint b) { return a = a * b; }

int N, Q;
vector<int> G[100000];
int K;

int H[100000];
int P[100000];
int V[100000];
int D[100000];
int S[100000];

void dfs(int u, int p) {
    S[u] = 1;
    for (int v : G[u]) if (v != p) {
        D[v] = D[u] + 1;
        dfs(v, u);
        S[u] += S[v];
    }
}

void dfs2(int u, int p, int h) {
    V[u] = K++;
    H[u] = h;
    P[u] = p;
    for (int v : G[u]) if (v != p && S[u] < S[v] * 2) dfs2(v, u, h);
    for (int v : G[u]) if (v != p && S[u] >= S[v] * 2) dfs2(v, u, v);
}

int lca(int u, int v) {
    for (;;) {
        if (V[u] > V[v]) swap(u, v);
        if (H[u] == H[v]) return u;
        v = P[H[v]];
    }
}

int BIT[100001];

void update(int k, int v) {
    for (int i = k + 1; i <= 100000; i += i & -i) {
        BIT[i] += v;
    }
}

int query(int k) {
    int res = 0;
    for (int i = k + 1; i > 0; i -= i & -i) {
        res += BIT[i];
    }
    return res;
}

int query(int l, int r) {
    return query(r) - query(l - 1);
}

int query_path(int u, int v) {
    int res = 0;
    for (;;) {
        if (V[u] > V[v]) swap(u, v);
        if (H[u] == H[v]) {
            res += query(V[u], V[v]);
            return res;
        } else {
            res += query(V[H[v]], V[v]);
            v = P[H[v]];
        }
    }
}

mint dp0[310];
mint dp1[310];

int main(void) {
    scanf("%d %d", &N, &Q);
    for (int i = 0; i < N - 1; i++) {
        int u, v;
        scanf("%d %d", &u, &v);
        u--;
        v--;
        G[u].push_back(v);
        G[v].push_back(u);
    }
    dfs(0, -1);
    dfs2(0, -1, 0);
    while (Q--) {
        int k, m, r;
        scanf("%d %d %d", &k, &m, &r);
        r--;
        vector<pair<int, int>> a(k);
        for (int i = 0; i < k; i++) {
            int x;
            scanf("%d", &x);
            x--;
            a[i].first = D[r] + D[x] - 2 * D[lca(r, x)];
            a[i].second = x;
        }
        sort(a.begin(), a.end());
        for (int i = 0; i <= m; i++) {
            dp0[i] = 0;
        }
        dp0[0] = 1;
        for (auto p : a) {
            for (int i = 0; i <= m; i++) {
                dp1[i] = dp0[i];
                dp0[i] = 0;
            }
            int cnt = query_path(r, p.second);
            for (int i = 0; i <= m; i++) {
                dp0[i + 1] += dp1[i];
                dp0[i] += dp1[i] * max(0, i - cnt);
            }
            update(V[p.second], 1);
        }
        for (auto p : a) {
            update(V[p.second], -1);
        }
        mint ans;
        for (int i = 0; i <= m; i++) {
            ans += dp0[i];
        }
        printf("%d\n", ans.n);
    }
    return 0;
}

Lunar New Year and a Recursive Sequence

https://codeforces.com/contest/1106/problem/F

リサンタイスウについて かく。ホウ 7 では 3 のベキをならべると $[1,3,3^2,3^3,3^4,3^5]=[1,3,2,6,4,5]$ になるけど 1,2,3,4,5,6 が すべて あらわれている。こういうかずは ゲンシコンだとかセイセイゲンだとか よばれている。ソスウをホウとしたとき ゲンシコンは かならずある。(これを しめすのは ちょっぴり むずかしい)

3 をなんジョウすると x になるのかを log x とかくことにする。たとえば log 3=1, log 2=2, log 5=5, log 1=0 である。このとき

\begin{align}
\log (ab) \equiv \log(a) + \log(b) \pmod{p - 1}
\end{align}

がなりたつ(isomorphism だね)。これによって このモンダイのゼンカシキは とてもカンタンになる。たとえば $a_{n+2}=a_{n+1}a_{n}$ みたいなゼンカシキを $\log a_{n+2} = \log a_{n+1} + \log a_{n}$ みたいにできる。のこりのパートは よくあるやつなので はぶきます。

リサンタイスウをもとめるアルゴリズムとしては baby step giant step があるけど セツメイは ほかにまかせます。

#include <stdio.h>
#include <vector>
#include <map>
#include <algorithm>

using namespace std;

#define MOD 998244353

long long mpow(long long a, long long b, long long mod) {
    long long res = 1;
    while (b > 0) {
        if (b & 1) res = res * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return res;
}

int mlog(long long r, long long b, long long mod) {
    map<long long, int> mp;
    for (int i = 0; i < 100000; i++) {
        mp[mpow(r, i, mod)] = i;
    }
    for (int i = 0; i < mod - 1; i += 100000) {
        if (mp.count(b * mpow(r, mod - 1 - i, mod) % mod)) {
            return i + mp[b * mpow(r, mod - 1 - i, mod) % mod];
        }
    }
    return -1;
}

vector<long long> pmod(vector<long long> a, vector<long long> p, long long mod) {
    const int n = a.size();
    const int m = p.size();
    for (int i = n - 1; i >= m; i--) {
        for (int j = 0; j < m; j++) {
            a[i - m + j] += a[i] * p[j];
            a[i - m + j] %= mod;
        }
    }
    a.resize(m);
    return a;
}

vector<long long> pmul(vector<long long> a, vector<long long> b, vector<long long> p, long long mod) {
    const int n = a.size();
    const int m = b.size();
    vector<long long> c(n + m - 1);
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < m; j++) {
            c[i + j] += a[i] * b[j];
            c[i + j] %= mod;
        }
    }
    return pmod(c, p, mod);
}

vector<long long> ppow(vector<long long> a, long long b, vector<long long> p, long long mod) {
    vector<long long> res(1, 1);
    while (b > 0) {
        if (b & 1) res = pmul(res, a, p, mod);
        a = pmul(a, a, p, mod);
        b >>= 1;
    }
    return res;
}

long long extgcd(long long a, long long b, long long c, long long mod) {
    long long s = 1;
    long long t = 0;
    while (b != 0) {
        long long q = a / b;
        a -= b * q;
        s -= t * q;
        swap(a, b);
        swap(s, t);
    }
    if (c % a != 0) return -1;
    if (s < 0) s += mod;
    return c / a * s % mod;
}

int main() {
    int K;
    long long N, M;
    scanf("%d", &K);
    vector<long long> p(K);
    for (int i = K - 1; i >= 0; i--) {
        scanf("%lld", &p[i]);
        p[i] %= MOD - 1;
    }
    scanf("%lld %lld", &N, &M);
    N--;
    M = mlog(3, M, MOD);
    vector<long long> a(2);
    a[1] = 1;
    a = ppow(a, N, p, MOD - 1);
    long long ans = extgcd(a[K - 1], MOD - 1, M, MOD - 1);
    if (ans == -1) {
        puts("-1");
    } else {
        ans = mpow(3, ans, MOD);
        printf("%lld\n", ans);
    }
    return 0;
}

そういえば カクチョウユークリッドのゴジョホウなんだけど たとえば $13x+8y=1$ を とくばあい

\begin{gather}
13\times 1 + 8 \times 0 = 13 \\
13 \times 0 + 8 \times 1 = 8 \\
13 \times 1 + 8 \times (-1) = 5 \\
13 \times (-1) + 8 \times 2 = 3 \\
13 \times 2 + 8 \times (-3) = 2 \\
13 \times (-3) + 8 \times 5 = 1 \\
13 \times 8 + 8 \times (-13) = 0
\end{gather}

みたいに できる。(L3)=(L1)-(L2), (L4)=(L3)-(L2), (L5)=(L4)-(L3) ... みたいな。むかしは もっとメンドウなこと かんがえてた。

$A^K - b_1 A^{K - 1} - \cdots - b_K E = O$ なら $A^n=A^n \bmod (A^K - b_1 A^{K - 1} - \cdots - b_K E)$ も まえにかいたキがするけど このカンケイをつかうと タコウシキに はなしが おきかえられる。