[jzoj 6084] [GDOI2019模拟2019.3.25] 礼物 [luogu 4916] 魔力环解题报告(莫比乌斯反演+生成函数)

题目链接：

https://jzoj.net/senior/#main/show/6084

https://www.luogu.org/problemnew/show/P4916

题目：

题解：

注：本题解大部分摘自Imagine大佬提供在洛谷的题解

我们设$f(x)$表示最小循环节长度为x的合法序列数，那么有$ans=\sum_{d|gcd(n,m)}\frac{1}{d}f(d)$

这是因为最小循环节为d的序列对应的环会被计算d次，比如

0101，最小循环节长度为 2（循环节为 01），其对应的环会被统计 2 次（序列 0101 与 1010）
0111，最小循环节长度为 4（循环节为 0111），其对应的环会被统计 4 次（序列 0111，1011，1101，1110）
1111，最小循环节长度为 1（循环节为 1），其对应的环会被统计 1 次（序列 1111）

直接求$f(x)$并不好求，我们令$g(x)=\sum_{d|x}f(d)$，即$g(x)$表示最小循环节是x的因子的合法序列的个数。那么可以莫比乌斯反演来求$f(x)$,$f(x)=\sum_{d|x}\mu(d)g(\frac{x}{d})$

问题转化为了如何快速的求$g(x)$

注意到这样一个事实：对于长度为 $a$ 且恰好有 $b$ 颗黑色珠子的一段珠子，我们求 $g(a)$，等同于求下面方程的整数解的数量：

$$x_0 + x_1 + ... + x_{a - b} = b(0 \leq x_i \leq k, 0 \leq x_0 + x_{a - b} \leq k)$$

即被 $a - b$ 颗白色珠子划分开的 $a - b + 1$ 段黑色珠子的和为 $b$，且满足每连续一段长度不超过 $k$ 的限制条件。运用生成函数的知识，求上面方程的解的数量等同于求如下多项式 $h(x)$ 中 $x^b$ 的系数：

$$h(x) = \left(\sum_{i = 0}^{k} x^i\right) ^ {a - b - 1} \left( \left(\sum_{i = 0}^{k} x^i\right)^2{\rm mod}\ x^{k + 1}\right)$$

进一步地，有

$$h(x) = \left(\sum_{i = 0}^{k} x^i\right) ^ {a - b - 1} \left(\sum_{i = 0}^{k} (i +1)x^i\right)$$

我们转化一下。由于 $\sum_{i = 0}^k x^i = \frac{1 - x^{k + 1}}{1 - x}$，因此有：

$$h(x) = \left(\frac{1 - x^{k + 1}}{1 - x}\right) ^ {a - b - 1} \left(\sum_{i = 0}^{k} (i +1)x^i\right)$$

再展开右侧的式子 $\sum_{i = 0}^k(i + 1)x^i$：

$$\begin{aligned}\sum_{i = 0}^k (i +1)x^i &= x^0 + 2x^1 + 3x^2 + \cdots + (k + 1)x^k\\ &= (x^0 + x^1 + \cdots + x^k) + (x^1 +x^2 + \cdots + x^k)+ \cdots + x^k \\ &= \frac{x^0 - x^{k + 1}}{1 - x} + \frac{x^1 - x^{k + 1}}{1 - x} + \cdots + \frac{x^k - x^{k - 1}}{1 - x} \\ &= \frac{(x^0 + x^1 + \cdots + x^k) - (k + 1)x^{k + 1}}{1 - x} \\ &= \frac{\frac{x^0 - x^{k + 1}}{1 - x} - (k + 1)x^{k + 1}}{1 - x} \\ &= \frac{1 - (k + 2)x^{k + 1} + (k + 1)x^{k + 2}}{(1 - x)^2}\end{aligned}$$

因此，我们得到了：

$$\begin{aligned}h(x) &= \left(\frac{1 - x^{k + 1}}{1 - x}\right) ^ {a - b - 1} \frac{1 - (k + 2)x^{k + 1} + (k + 1)x^{k + 2}}{(1 - x)^2} \\ &= \frac{(1 - x^{k + 1})^{a - b - 1}}{(1 - x)^{a - b + 1}}(1 - (k + 2)x^{k + 1} + (k + 1)x^{k + 2})\end{aligned}$$

其中，$(1 - x^{k + 1})^{a - b - 1}$ 可化为 $\sum_{i = 0}^{\infty}\binom{a - b - 1}{i}(-1)^ix^{(k + 1)i}$，而 $\frac{1}{(1 - x)^{a - b + 1}}$ 即 $(1 - x)^{-(a - b + 1)}$，可通过负整数次幂的二项式定理化为 $\sum_{i = 0}^{\infty}\binom{a - b + i}{i}x^i$，因此有：

$$h(x) = \left(\sum_{i = 0}^{\infty}\binom{a - b - 1}{i}(-1)^ix^{(k + 1)i}\right)\left(\sum_{i = 0}^{\infty}\binom{a - b + i}{i}x^i\right)(1 - (k + 2)x^{k + 1} + (k + 1)x^{k + 2})$$

当把 $h(x)$ 化成该形式后，要求 $h(x)$ 中 $x^b$ 的系数就变得非常简单了。记 $s_1 = \sum_{(k + 1)i + j = b}(-1)^i\binom{a - b - 1}{i}\binom{a - b+ j}{j}$，$s_2 = (k + 2)\sum_{(k + 1)i + j = b - k - 1}(-1)^i\binom{a - b - 1}{i}\binom{a - b+ j}{j}$，$s_3 = (k + 1)\sum_{(k + 1)i + j = b - k - 2}(-1)^i\binom{a - b - 1}{i}\binom{a - b+ j}{j}$，$x^b$ 的系数为 $w$，那么有 $w = s_1 - s_2 + s_3$。

求 $s_1, s_2, s_3$ 只需按照 $s_1, s_2, s_3$ 的式子枚举 $i$ 即可，因为 $i$ 确定 $j$ 也就确定了。因此，我们可以在 $\frac{b}{k + 1}$ 的时间内求出 $h(x)$ 中 $x^b$ 的系数。

除去反演部分，我们就能够在 $\frac{\sigma(n)}{k + 1}$ 的时间内解决此题，其中，$\sigma(n)$ 表示 $n$ 的约数和。由于 $\sigma(n)$ 可近似看作 $n\ {\rm log}\ {\rm log}\ n$，接近线性，因此时间复杂度是非常优秀的。

#include<algorithm>

#include<cstring>

#include<cstdio>

#include<iostream>

using namespace std;

const int N=1e6+;

const int mo=;

int k,cnt;

int mu[N],prime[N],vis[N];

int f[N],g[N],fac[N],inv[N];

inline int read()

{

    char ch=getchar();int s=,f=;

    while (ch<''||ch>'') {if (ch=='-') f=-;ch=getchar();}

    while (ch>=''&&ch<='') {s=(s<<)+(s<<)+ch-'';ch=getchar();}

    return s*f;

}

int qpow(int a,int b)

{

    int re=;

    for (;b;b>>=,a=1ll*a*a%mo) if (b&) re=1ll*re*a%mo;

    return re;

}

void pre()

{

    mu[]=;

    for (int i=;i<N;i++)

    {

        if (!vis[i])

        {

            prime[++cnt]=i;

            mu[i]=-;

        }

        for (int j=;j<=cnt&&prime[j]*i<N;j++)

        {

            vis[prime[j]*i]=;

            if (i%prime[j]) mu[prime[j]*i]=-mu[i];

            else break;

        }

    }

    fac[]=;

    for (int i=;i<N;i++) fac[i]=1ll*fac[i-]*i%mo;

    inv[N-]=qpow(fac[N-],mo-);

    for (int i=N-;i>=;i--) inv[i]=1ll*inv[i+]*(i+)%mo;

}

int gcd(int a,int b)

{

    if (!b) return a;

    return gcd(b,a%b);

}

int C(int a,int b)

{

    return 1ll*fac[a]*inv[b]%mo*inv[a-b]%mo;

}

void add(int &a,int b)

{

    if (b<) b+=mo;

    a=(a+b)%mo;

}

int calc(int n,int m)

{

    int res=;

    for (int i=;i*(k+)<=m;++i)

    {

        int j=m-i*(k+);

        if (i&) add(res,1ll*C(n-m-,i)*C(n-m+j,j)%mo*(mo-)%mo);

        else add(res,1ll*C(n-m-,i)*C(n-m+j,j)%mo*%mo);

        j=m-i*(k+)-k-;

        if (j>=)

        {

            if (i&) add(res,1ll*(k+)*C(n-m-,i)%mo*C(n-m+j,j)%mo*%mo);

            else add(res,1ll*(k+)*C(n-m-,i)%mo*C(n-m+j,j)%mo*(mo-)%mo);

        }

        j=m-i*(k+)-k-;

        if (j>=)

        {

            if (i&) add(res,1ll*(k+)*C(n-m-,i)%mo*C(n-m+j,j)%mo*(mo-)%mo);

            else add(res,1ll*(k+)*C(n-m-,i)%mo*C(n-m+j,j)%mo*%mo);

        }

    }

    return res;

}

int n,m;

int main()

{

    freopen("gift.in","r",stdin);

    freopen("gift.out","w",stdout);

    pre();

    int T=read();

    while (T--)

    {

        memset(f,,sizeof(f));

        memset(g,,sizeof(g));

        n=read();m=read();k=read();

        int d=gcd(n,m);

        for (int i=;i<=d;i++)

        {

            if (d%i) continue;

            g[n/i]=calc(n/i,m/i);

        }

        for (int i=;i<=n;i++)

            for (int j=i;j<=n;j+=i)

                add(f[j],mu[i]*g[j/i]);

        int ans=;

        for (int i=;i<=n;i++)

        {

            if (n%i) continue;

            add(ans,1ll*f[i]*qpow(i,mo-)%mo);

        }

        printf("%d\n",ans);

    }

    return ;

}

巴特西

[jzoj 6084] [GDOI2019模拟2019.3.25] 礼物 [luogu 4916] 魔力环解题报告(莫比乌斯反演+生成函数)

题目链接：

题目：

题解：

最新文章

热门文章

巴特西

[jzoj 6084] [GDOI2019模拟2019.3.25] 礼物 [luogu 4916] 魔力环 解题报告(莫比乌斯反演+生成函数)

题目链接：

题目：

题解：

最新文章

热门文章

[jzoj 6084] [GDOI2019模拟2019.3.25] 礼物 [luogu 4916] 魔力环解题报告(莫比乌斯反演+生成函数)