Project Euler 560

题目

Coprime Nim

Coprime Nim is just like ordinary normal play Nim, but the players may only remove a number of stones from a pile that is coprime with the current size of the pile. Two players remove stones in turn. The player who removes the last stone wins.

Let \(L(n, k)\) be the number of losing starting positions for the first player, assuming perfect play, when the game is played with \(k\) piles, each having between \(1\) and \(n - 1\) stones inclusively.

For example, \(L(5, 2) = 6\) since the losing initial positions are \((1, 1), (2, 2), (2, 4), (3, 3), (4, 2)\) and \((4, 4)\).

You are also given \(L(10, 5) = 9964, L(10, 10) = 472400303, L(10^3, 10^3)  \bmod 1 000 000 007 = 954021836\).

Find \(L(10^7, 10^7) \bmod 1 000 000 007\)

快速沃尔什变换

目标：给定两个长度为\(N=2^n\)的序列\(a,b\)（两个序列从下标从\(0\)开始），求一个长度为\(N\)的序列\(c\)，满足：

\[c_k=\sum_{i \circ j=k} a_i\cdot b_j\]

其中运算符\(\circ\)是三种位运算的运算符之一：与运算\(\&\)，或运算\(|\)，异或运算\(\oplus\)。

快速沃尔什变换可以将整个计算过程的时间复杂度从\(O(N^2)\)下降到\(O(N \log N)\)。

此处不详述这个算法，此处为参考文章和代码。

解决方案

令\(N=10^7,K=10^7\)。

可以看出这是一个公平组合游戏，我们将每一堆石头视为独立的一个游戏局面，使用SG定理来解决本题。这里的石头一共有\(K\)堆，第\(i\)堆一共有\(a_i\)个，那么整个游戏的\(sg\)函数值为\(sg(a_1)\oplus sg(a_2)\oplus\dots sg(a_k)\)。

考虑目前的某一堆石头，一共有\(n\)个，每一次操作是取出\(d\)个石头，其中\(\gcd(n,d)=1\)。因此，可以写出一堆石头的\(sg\)函数定义：

\[ sg(n)= \left \{\begin{aligned} &0 & & \text{if}\quad n=0 \\ &\text{mex}(\{sg(n-d)|\gcd(i,d)=1\}) & & \text{else} \end{aligned}\right. \]

其中运算\(\text{mex}(s)\)表示集合\(s\)中未出现的最小的非负整数。

那么考虑运行这一份代码，打印出前一部分项的\(sg\)函数。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=100,M=100;
int sg[N+4];
bool mex[M+4];
int pr[N/10+100],v[N+4];
int main(){
    for(int i=1;i<=N;i++){
        memset(mex,0,sizeof(mex));
        for(int j=1;j<=i;j++)
            if(__gcd(i,j)==1) mex[sg[i-j]]=1;
        int j=0;
        for(;mex[j];++j);
        sg[i]=j;
        printf("%d %d\n",i,sg[i]);
    }
}

发现规律如下：

• \(sg(1)=1\)
• 如果\(n\)为偶数，那么\(sg(n)=0\)
• 否则，假设\(n\)的最小质因子为\(p\)，如果质数表中第\(k\)个质数（\(k\)从\(1\)开始算起）为\(p\)，那么\(sg(n)=k\)。

最终问题转化为有多少种方法可以将\(sg\)函数值填入一个长度为\(K\)的序列，并且这\(K\)个值的异或和为\(0\)。

令状态\(g(i,j)(0\le i\le K,j\ge 0)\)表示填入长度为\(i\)的序列，并且异或和为\(j\)的方案数有多少个。那么不难写出状态转移方程为：

\[ g(i,j)= \left \{\begin{aligned} &1 & & \text{if}\quad i=0\land j=0 \\ &0 & & \text{else if}\quad i=0 \\ &\sum_{k} g(i-1,j\oplus k) \cdot c[k] & & \text{else} \end{aligned}\right. \]

注意到，此时\(sg\)函数值比较大，可以达到\(O\left(\dfrac{N}{\log N}\right)\)的级别。对于\(Q\)这种数据范围而言，直接进行转移效率非常低。

如果我们已经求出了长度为\(a\)的序列的各种填法\(g(a,\cdot)\)和长度为\(b\)的各种填法\(g(b,\cdot)\)，不难知道我们可以通过卷积组合出\(g(a+b,\cdot)\)。

枚举所有的\(g(a,i)\)和\(g(b,j)\)，将这两部分值直接合并起来：

\(g(a,i)\cdot g(b,j)\rightarrow g(a+b,i\oplus j)\)

直接枚举这两部分的值效率非常低，可以使用快速沃尔什变换加速计算这两部分的合并。

因此，这给了我们一个方案：依次求出\(g(2^0,\cdot),g(2^1,\cdot),g(2^2,\cdot),\dots\)。然后针对\(Q\)，选择这些求出的\(g(2^i,\cdot)\)进行合并即可。

最终答案为\(g(Q,0)\)。

代码

为了加速，特地设计了一个powXor的函数。如果直接按照上面的快速沃尔什变换模板使用，那么总共需要做\(O(\log Q)\)次变换和逆变换。而此处由于是多次求值，因此总共只需要\(3\)次变换和逆变换。

# include <bits/stdc++.h>
# include "tools.h"
using namespace std;
typedef long long ll;
const int N=10000000,M=10000000;
int mod= 1000000007;
int pr[N/5+100],v[N+4],m=0;
int main(){
    for(int i=2;i<=N;i++){
        if(v[i]==0){
            pr[++m]=i;v[i]=m;
        }
        for(int j=1;j<=m;j++){
            if(j>v[i]||pr[j]>N/i) break;
            v[i*pr[j]]=j;
        }
    }
    v[1]=1;
    for(int i=2;i<=N;i+=2)
        v[i]=0;
    vector<int>b(m+1);
    for(int i=1;i<N;i++)
        ++b[v[i]];
    fwt.init_mod(mod);
    int ans=fwt.powXor(b,M)[0];
    printf("%d\n",ans);
}