51nod 1778 小Q的集合 lucas定理+线性筛+数学

题意

小Q有一个集合$S$，它的元素个数$|S|=n$。
对于$S$的任意一个子集合$T$，定义$f(T)=|T|^k$，定义$T$关于$S$的补集为$S−T$。
小Q想知道，如果他等概率地选择一个$S$的子集$T$，那么$f(T)−f(S−T)$的方差是多少。
由于这个方差值可能很大，不妨设其为$v$，你只需要给出$(v⋅2^n)\bmod m$的值即可。
$k ≤ n ≤ 10^{10^6},1<=k<=10^6,m是质数,2<=m<=10^6$

分析

注意到$f(T)-f(S-T)$的平均数是0，那么我们要求的其实就是

$$\sum_{i=0}^nC_n^i(i^k-(n-i)^k)^2$$
因为根据lucas定理有 $$C_n^i\bmod m=C_{n/m}^{i/m}*C_{n\bmod m}^{i\bmod m}\bmod m$$ 所以我们只用考虑所有满足 $i\bmod m<=n\bmod m$的 $i$即可。
然后就开始推式子：
$$ans=\sum_{i=0}^{\lfloor\frac{n}{m}\rfloor}\sum_{j=0}^{n\bmod m}C_n^{im+j}((im+j)^k-(n-im-j)^k)^2$$
$$=\sum_{i=0}^{\lfloor\frac{n}{m}\rfloor}C_{\lfloor\frac{n}{m}\rfloor}^i\sum_{j=0}^{n\bmod m}C_{n\bmod m}^j(j^k-(n-j)^k)^2$$
$$=2^{\lfloor\frac{n}{m}\rfloor}\sum_{j=0}^{n\bmod m}C_{n\bmod m}^j(j^k-(n-j)^k)^2$$
前面2的指数可以用费马小定理来优化，后面可以线性筛处理自然数幂，线性预处理组合数，这样总的复杂度就是 $O(n)$的了。

代码

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
using namespace std;

typedef long long LL;

const int N=1000005;
const int maxn=1000000;

int m,k,s[N],prime[N],tot,a[N],jc[N],ny[N];
bool not_prime[N];
char str[N];

int ksm(int x,int y)
{
    int ans=1;
    while (y)
    {
        if (y&1) ans=(LL)ans*x%m;
        x=(LL)x*x%m;y>>=1;
    }
    return ans;
}

int C(int x,int y)
{
    return (LL)jc[x]*ny[y]%m*ny[x-y]%m;
}

int sqr(int x)
{
    return (LL)x*x%m;
}

int divi(int c)
{
    int s,g=0;
    for (int i=1;i<=maxn;i++)
    {
        s=g*10+a[i];
        a[i]=s/c;
        g=s%c;
    }
    return g;
}

void get_prime(int n)
{
    s[1]=1;
    for (int i=2;i<=n;i++)
    {
        if (!not_prime[i]) prime[++tot]=i,s[i]=ksm(i,k);
        for (int j=1;j<=tot&&i*prime[j]<=n;j++)
        {
            not_prime[i*prime[j]]=1;
            s[i*prime[j]]=(LL)s[i]*s[prime[j]]%m;
            if (i%prime[j]==0) break;
        }
    }
}

int main()
{
    scanf("%s%d%d",str+1,&k,&m);int len=strlen(str+1);
    for (int i=1;i<=len;i++) a[maxn-len+i]=str[i]-'0';
    int r1=divi(m),r2=divi(m-1);
    get_prime(r1);
    jc[0]=jc[1]=ny[0]=ny[1]=1;
    for (int i=2;i<=r1;i++) jc[i]=(LL)jc[i-1]*i%m,ny[i]=(LL)(m-m/i)*ny[m%i]%m;
    for (int i=2;i<=r1;i++) ny[i]=(LL)ny[i-1]*ny[i]%m;
    int ans=0;
    for (int i=0;i<=r1;i++) (ans+=(LL)C(r1,i)*sqr(s[i]-s[r1-i])%m)%=m;
    ans=(LL)ans*ksm(2,r2)%m;
    printf("%d",ans);
    return 0;
}