【BZOJ4589】Hard Nim-快速幂+FWT

测试地址：Hard Nim
做法：本题需要用到快速幂+FWT。
我们知道Nim游戏先手必败的条件为，所有堆中的石子数异或和为0$0$，于是我们就是要求石子数异或和为0$0$的方案数。我们直观上感觉这个东西很像生成函数的卷积，但是又不是一般的卷积，因为这里并不是石子数的和，而是异或和。
这时我们要重新定义卷积。令向量c$c$为向量a$a$和b$b$关于运算⊕$\oplus$的卷积，则c$c$为：
ck=∑i⊕j=kai⋅bj$c_k=\sum _{i\oplus j=k}{a}_i\cdot b_j$
根据这个定义，我们平时所说的卷积实际上应该是关于加法的卷积，而在这里我们就要求关于异或的卷积。于是我们就引出了FWT，即快速沃尔什变换。
快速沃尔什变换提供了O(nlogn)$O(n\log n)$（其中n$n$为向量维数，也即长度）求出位运算卷积的方法，包含关于与、或、异或三种位运算的卷积的计算，不同运算有不同公式，但本质上类似。具体的对FWT的理解，因为网上有其他大佬的理解，这里我先不写，等有时间再进行总结。
于是回到这题，实际上就是求n$n$个相同向量的异或卷积，用快速幂和FWT即可做到O(mlogmlogn)$O(m\log m\log n)$。
以下是本人代码：

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod=1000000007;
const ll inv=500000004;
int n,m,prime[70010];
ll s[70010],now[70010];
bool vis[70010];

void calc_prime()
{
    prime[0]=0;
    memset(vis,0,sizeof(vis));
    for(int i=2;i<=m;i++)
    {
        if (!vis[i]) prime[++prime[0]]=i;
        for(int j=1;j<=prime[0]&&i*prime[j]<=m;j++)
        {
            vis[i*prime[j]]=1;
            if (i%prime[j]==0) break;
        }
    }
}

void FWT(ll *a,int type)
{
    for(int mid=1;mid<m;mid<<=1)
        for(int l=0;l<m;l+=(mid<<1))
            for(int k=0;k<mid;k++)
            {
                ll x=a[l+k],y=a[l+mid+k];
                a[l+k]=(x+y)%mod;
                a[l+mid+k]=(x-y+mod)%mod;
                if (type==-1)
                {
                    a[l+k]=a[l+k]*inv%mod;
                    a[l+mid+k]=a[l+mid+k]*inv%mod;
                }
            }
}

void power(int n)
{
    FWT(s,1),FWT(now,1);
    while(n)
    {
        if (n&1)
        {
            for(int i=0;i<=m;i++)
                s[i]=s[i]*now[i]%mod;
        }
        for(int i=0;i<=m;i++)
            now[i]=now[i]*now[i]%mod;
        n>>=1;
    }
    FWT(s,-1);
}

int main()
{
    m=50000;
    calc_prime();

    while(scanf("%d%d",&n,&m)!=EOF)
    {
        memset(now,0,sizeof(now));
        memset(s,0,sizeof(s));
        s[0]=1;
        for(int i=1;i<=prime[0]&&prime[i]<=m;i++)
            now[prime[i]]=1;

        int x=1;
        while(x<=m) x<<=1;
        m=x;
        power(n);

        printf("%lld\n",s[0]);
    }

    return 0;
}