Solution 「NOI Simu.」树

\(\mathscr{Description}\)

  给定 \(n\), 集合 \(\{a_m\}\), 称一棵无标号但儿子有序的有根树合法, 当且仅当叶子权值存在一个因数在 \(\{a_m\}\) 中, 且叶子点权积不超过 \(n\). 此外, 要求不存在仅含一个儿子的结点. 求有多少棵合法树. 答案模 \(10^9+7\).

\(\mathscr{Solution}\)

  令 \(L(x)=\sum_{i}[\exists t\in\{a_m\},~t\mid i]x^i\), \(F(x)\) 为有根树关于叶子乘积和的 GF, 那么

\[ F(x)=L(x)+\sum_{i\ge 2}F^i(x), \]

其中卷积是 Dirichlet 卷积. 推一下式子:

\[ F(x)-\frac{F^2(x)}{\epsilon-F(x)}=L(x)\\ \Rightarrow 2F^2(x)-(L(x)+\epsilon)F(x)+L(x)=0\\ \Rightarrow F(x)=\frac{L(x)+\epsilon-\sqrt{L^2(x)-6L(x)+\epsilon}}{4}. \]

  我们想要求 \(F(x)\) 前缀系数和, 如果能求出 \(L^2(x)-6L(x)+\epsilon\) 的 \(\sqrt n\) 个前缀和, 并以此算出开根后 \(n\) 处的前缀和, 实际上就结束了.

  不太平凡的只有 Dirichlet 开根, 令 \(G^2(x)=T(x)=L^2(x)-6L(x)+\epsilon\), 那么

\[ \sum_{i=1}^m\sum_{d\mid i}g_dg_{i/d}=\sum_{i=1}^mt_i\\ \Rightarrow 2S_g(m)=S_t(m)-\sum_{i=2}^{m/2}g_iS_g(m/i)+S_g(m/2+1). \]

  经过暴力卷积预处理分块, 可以做到 \(\mathcal O(n^{2/3}\ln^{1/3}n)\).

\(\mathscr{Code}\)

/*+Rainybunny+*/

#include <bits/extc++.h>

#define rep(i, l, r) for (int i = l, rep##i = r; i <= rep##i; ++i)
#define per(i, r, l) for (int i = r, per##i = l; i >= per##i; --i)

typedef long long LL;

const int MOD = 1e9 + 7, MAXM = 8, INV2 = 500000004, INV4 = 250000002;
const int THRES = 2e6, MAXSQRT = 1e5; // (n/ln(n))^{2/3}.
// const int THRES = 5;
LL n, lcm[1 << MAXM];
int sn, m, a[MAXM + 5], L[THRES + 5], L2[THRES + 5], G[THRES + 5];
int pn, pr[THRES + 5], mp[THRES + 5];
bool npr[THRES + 5];

inline int mul(const int u, const int v) { return 1ll * u * v % MOD; }
inline void subeq(int& u, const int v) { (u -= v) < 0 && (u += MOD); }
inline int sub(int u, const int v) { return (u -= v) < 0 ? u + MOD : u; }
inline void addeq(int& u, const int v) { (u += v) >= MOD && (u -= MOD); }
inline int add(int u, const int v) { return (u += v) < MOD ? u : u - MOD; }
inline int mpow(int u, int v) {
    int ret = 1;
    for (; v; u = mul(u, u), v >>= 1) ret = mul(ret, v & 1 ? u : 1);
    return ret;
}

inline void init() {
    rep (i, 2, THRES) {
        if (!npr[i]) mp[pr[++pn] = i] = i;
        for (int j = 1, t; j <= pn && (t = i * pr[j]) <= THRES; ++j) {
            npr[t] = true, mp[t] = pr[j];
            if (!(i % pr[j])) break;
        }
    }

    rep (S, 1, (1 << m) - 1) {
        lcm[S] = 1;
        rep (i, 0, m - 1) if (S >> i & 1) {
            lcm[S] = lcm[S] / std::__gcd(lcm[S], 0ll + a[i]) * a[i];
        }
        if (__builtin_parity(S)) rep (i, 1, THRES / lcm[S]) ++L[lcm[S] * i];
        else rep (i, 1, THRES / lcm[S]) --L[lcm[S] * i];
    }
    rep (i, 1, THRES) rep (j, 1, THRES / i) addeq(L2[i * j], mul(L[i], L[j]));

    G[1] = 1;
    rep (i, 2, THRES) {
        static int pc, pv[30], alp[30]; pc = 0;
        for (int t = i; t > 1;) {
            pv[++pc] = mp[t], alp[pc] = 0;
            while (mp[t] == pv[pc]) ++alp[pc], t /= pv[pc];
        }

        std::function<int(int, int)>
        contri = [&](const int x, int v)->int {
            if (x > pc) return mul(G[v], G[i / v]);
            int ret = 0;
            rep (j, 1, alp[x]) addeq(ret, contri(x + 1, v)), v *= pv[x];
            return add(ret, contri(x + 1, v));
        };
        G[i] = mul(sub(sub(L2[i], mul(6, L[i])), contri(1, 1)), INV2);
    }

    rep (i, 1, THRES) {
        addeq(L[i], L[i - 1]), addeq(L2[i], L2[i - 1]), addeq(G[i], G[i - 1]);
    }
}

inline int id(const LL x) { return x <= sn ? x : sn + n / x; }

inline int calcL(const LL m) {
    static int mem[MAXSQRT * 2 + 5];
    static int vis[MAXSQRT * 2 + 5];
    if (m <= THRES) return L[m];
    int h = id(m);
    if (vis[h]) return mem[h];
    int& ret = mem[h]; vis[h] = true;
    rep (S, 1, (1 << ::m) - 1) {
        (__builtin_parity(S) ? addeq : subeq)(ret, m / lcm[S] % MOD);
    }
    return ret;
}

inline int calcL2(const LL m) {
    static int mem[MAXSQRT * 2 + 5];
    static int vis[MAXSQRT * 2 + 5];
    if (m <= THRES) return L2[m];
    int h = id(m);
    if (vis[h]) return mem[h];
    int& ret = mem[h] = 0; vis[h] = true;
    for (LL l = 2, r; l <= m; l = r + 1) {
        r = m / (m / l);
        addeq(ret, mul(calcL(m / l), sub(calcL(r), calcL(l - 1))));
    }
    return ret;
}

inline int calcG(const LL m) {
    static int mem[MAXSQRT * 2 + 5];
    static int vis[MAXSQRT * 2 + 5];
    if (m <= THRES) return G[m];
    int h = id(m);
    if (vis[h]) return mem[h];
    int& ret = mem[h] = add(sub(calcL2(m), mul(6, calcL(m))), 1); vis[h] = 1;
    for (LL l = 2, r; l <= m; l = r + 1) {
        r = m / (m / l);
        if (m / l == 1) addeq(ret, calcG(l - 1));
        else subeq(ret, mul(calcG(m / l), sub(calcG(r), calcG(l - 1))));
    }
    ret = mul(ret, INV2);
    return ret;
}

int main() {
    freopen("tree.in", "r", stdin);
    freopen("tree.out", "w", stdout);

    scanf("%lld %d", &n, &m), sn = sqrt(n);
    rep (i, 0, m - 1) scanf("%d", &a[i]);

    init();

    int ans = mul(sub(add(calcL(n), 1), calcG(n)), INV4);
    printf("%d\n", ans);
    return 0;
}