「2018 集训队互测 Day 1」完美的集合 题解

「2018 集训队互测 Day 1」完美的集合 题解

对于选出的集合\(S\),满足条件的\(x\)一定构成一个联通块。

所以可以容斥，对于每一个点算出符合要求的集合数量，减去每一条边符合要求的集合数量。最终算的就是\({cnt\choose k}\mod 5^{23}\)。

类似\(ex\_lucas\)。

问题是算:

\[\frac{A!}{B!C!}\pmod {MOD} \]

由于会有5的倍数，所以需要将\(A!\)先除去所有5的因子。令\(f(A)=\frac{A!}{5^{\sum_{j\geq 1} \lfloor\frac{A}{5^j}\rfloor}}\)。

然后可以变成算\(\prod_{1\leq i\leq A,i\mod 5\neq 0}i\)。

这个东西可以想到分治计算，令\(f_A(x)=\prod_{1\leq i\leq A,i\mod 5\neq 0}(x+i)\)。

如果\(A=10z\)，可以将\(f_A\)分成\(f_{A/2}(x),f_{A/2}(x+A/2)\)。

最后我们只需要知道\(f_A(0)\),也就是要知道\(f_{A/2}(0),f_{A/2}(A/2)\)。如果\(A/2\mod 5=0\),则我们可以只保留多项式\(\mod x^{23}\).

/**
 *    author:  gary
 *    created: 02.01.2022 10:39:00
**/
#include <bits/stdc++.h>
#define rb(a,b,c) for(int a=b;a<=c;++a)
#define rl(a,b,c) for(int a=b;a>=c;--a)
#define rep(a,b) for(int a=0;a<b;++a)
#define LL long long
#define PB push_back
#define POB pop_back
#define II(a,b) make_pair(a,b)
#define FIR first
#define SEC second
#define SRAND mt19937 rng(chrono::steady_clock::now().time_since_epoch().count())
#define random(a) rng()%a
#define ALL(a) a.begin(),a.end()
#define check_min(a,b) a=min(a,b)
#define check_max(a,b) a=max(a,b)
using namespace std;
const LL MOD = 11920928955078125;
const LL PHI = MOD / 5 * 4;
const int MAXN = 61;
LL quick(LL k1, LL k2) {
    LL k3 = 1;

    for (; k2; k2 >>= 1, k1 = ((__int128_t)(k1)) * k1 % MOD)
        if (k2 & 1)
            k3 = ((__int128_t)(k3)) * k1 % MOD;

    return k3;
}
LL inv(LL A) {
    return quick(A, PHI - 1);
}
LL count(LL x) {
    // x! 中5的个数
    LL cnt = 0;

    for (LL i = 5; i <= x; i *= 5)
        cnt += x / i;

    return cnt;
}
LL comb[1001][1001];
LL mul(LL A, LL B) {
    return ((__int128_t)(A)) * B % MOD;
}
vector<LL> mv(vector<LL> A, LL d) {
    // 向左平移d
    vector<LL> ans(23, 0);
    rep(i, 23) {
        rep(j, i + 1) {
            ans[j] += mul(mul(A[i], quick(d, i - j)), comb[i][j]);
            ans[j] %= MOD;
        }
    }
    return ans;
}
vector<LL> conv(vector<LL> A, vector<LL> B) {
    vector<LL> C(23, 0);
    rep(i, 23) rep(j, 23 - i)(C[i + j] += mul(A[i], B[j])) %= MOD;
    return C;
}
vector<LL> f(LL x) {
    // \prod_{1<=i<=x} (x+[i%5!=0]*i) mod x^23
    vector<LL> ans(23, 0);

    if (x == 0) {
        ans[0] = 1;
        return ans;
    }

    if (x % 10 != 0) {
        auto old = f(x - 1);
        LL coef = x;

        if (x % 5 != 0)
            rep(i, 23) {
            if (i != 22)
                (ans[i + 1] += old[i]) %= MOD;

            (ans[i] += mul(old[i], coef)) %= MOD;
        } else
            return old;

        return ans;
    }

    auto old = f(x / 2);
    return conv(old, mv(old, x / 2));
}
map<LL, LL> M;
LL calc(LL x) {
    // (x)!/(5^(count(x)))
    if (x < 5) {
        LL ans = 1;
        rb(i, 1, x) ans = ans * i;
        return ans;
    }

    if (M.find(x) != M.end())
        return M[x];

    __int128_t ans = calc(x / 5);
    M[x] = ans * (f(x)[0]) % MOD;
    return ans * (f(x)[0]) % MOD;
}
LL c(LL A, LL B) {
    if (A < B)
        return 0;

    // return comb[A][B];
    LL C = A - B;
    // A! / (B!*C!)
    LL cnt = count(A) - count(B) - count(C);

    if (cnt >= 23)
        return 0;

    return mul(mul(calc(A), inv(calc(B))), mul(inv(calc(C)), quick(5, cnt)));
}
int n, m, k;
LL Max;
int w[MAXN], v[MAXN];
vector<pair<int, int>> g[MAXN];
int dist[MAXN][MAXN];
bool ok[MAXN];
vector<int> dfn;
int siz[MAXN];
pair<LL, LL> dp[MAXN][10000 + 1];
void dfs(int now, int pre) {
    dfn.PB(now);
    siz[now] = 1;

    for (auto it : g[now])
        if (it.FIR != pre)
            dfs(it.FIR, now), siz[now] += siz[it.FIR];
}
void upd(pair<LL, LL> &A, pair<LL, LL> B) {
    if (B.FIR < A.FIR)
        return ;

    if (B.FIR > A.FIR)
        A = B;
    else
        A.SEC += B.SEC;
}
LL mx = -1e18;
LL solve(int x, int y) {
    dfn.clear();
    dfs(x, 0);

    if (w[x] > m)
        return 0;

    if (!ok[x] || !ok[y])
        return 0;

    rep(i, MAXN) rb(j, 0, 10000) dp[i][j] = II(-1e18, 0);
    dp[1][w[x]] = II(v[x], 1);
    rb(i, 1, n - 1) {
        rb(j, 0, m) {
            if (dp[i][j].SEC) {
                // cout<<i<<" "<<j<<" "<<dp[i][j].FIR<<" "<<dp[i][j].SEC<<" "<<dfn[i]<<" "<<y<<" "<<siz[i]<<endl;
                {
                    //choose
                    if (j + w[dfn[i]] <= m && ok[dfn[i]])
                        upd(dp[i + 1][j + w[dfn[i]]], II(dp[i][j].FIR + v[dfn[i]], dp[i][j].SEC));
                }
                {
                    //don't choose
                    if (dfn[i] != y) {
                        upd(dp[i + siz[dfn[i]]][j], dp[i][j]);
                    }
                }
            }
        }
    }
    pair<LL, LL> tmp = II(-1e18, 0);
    rb(j, 0, m) upd(tmp, dp[n][j]);
    check_max(mx, tmp.FIR);

    if (tmp.FIR != mx)
        return 0;

    return c(tmp.SEC, k);
}
int main() {
    ios::sync_with_stdio(false);
    cin.tie(0);
    rb(i, 0, 1000) comb[i][0] = 1;
    rb(i, 1, 1000) rb(j, 1, 1000) comb[i][j] = (comb[i - 1][j] + comb[i - 1][j - 1]) % MOD;
    cin >> n >> m >> k >> Max;
    rb(i, 1, n) cin >> w[i];
    rb(i, 1, n) cin >> v[i];
    rb(i, 1, n) rb(j, 1, n) dist[i][j] = 1e9;
    vector<pair<int, int>> edge;
    rb(i, 1, n - 1) {
        int A, B, C;
        cin >> A >> B >> C;
        edge.PB(II(A, B));
        g[A].PB(II(B, C));
        g[B].PB(II(A, C));
        dist[A][B] = dist[B][A] = C;
    }
    rb(i, 1, n) {
        rb(j, 1, n) ok[j] = 1;
        solve(i, i);
    }
    rb(i, 1, n) dist[i][i] = 0;
    rb(k, 1, n) rb(i, 1, n) rb(j, 1, n) check_min(dist[i][j], dist[i][k] + dist[k][j]);
    LL answer = 0;
    rb(i, 1, n) {
        memset(ok, 0, sizeof(ok));
        rb(j, 1, n) {
            ok[j] = (1ll * dist[i][j] * v[j] <= Max);
        }
        answer += solve(i, i);
        // cout<<answer<<endl;
        // return 0;
        answer %= MOD;
    }

    for (auto it : edge) {
        int x, y;
        tie(x, y) = it;
        memset(ok, 0, sizeof(ok));
        rb(j, 1, n) {
            ok[j] = (1ll * dist[x][j] * v[j] <= Max && 1ll * dist[y][j] * v[j] <= Max);
        }
        answer += MOD - solve(x, y);
        answer %= MOD;
    }

    cout << answer << endl;
    return 0;
}