题解 【SP422】

· · 题解

对于一个 1\sim n 的序列 a,如果想变成一个排列 p,我们连边 i\to p_i,显然会形成若干环,每个环的操作次数是 size-1。设环数为 x,那么答案就是 \sum_{i=1}^x size_i-1,而 \sum size_i=n,所以答案也可以写成 n-x。我们现在只需要求环数。

观察起始位置和中止位置的变化。先把所有下标变成 0-base 的,(i,j) 的初始位置是 2^bi+j,中止位置是 2^aj+i。容易发现,每个数的位置视为一个 a+b 位二进制数之后,并把这 a+b 个二进制位视为一个环,起始位置和中止位置的变化相当于把这个环旋转了 b 个单位(也可以说是旋转了 a 个单位,没有区别)。

如果 x 经过这个操作之后变成 y,那么我们会连边 x\to y。我们现在要数这个图上的环的数量。我们可以把 x_1\to x_2\to \cdots \to x_n\to x_1 这个环中的元素全部视为一个等价类,我们要数的就是等价类的数量。此时两个元素等价当且仅当他们可以通过 i(i\ge 0) 次「旋转 b 个单位」的操作变成相同的。

「旋转 b 个单位」的操作会把 a+b 个二进制位划分成 \gcd(a,b) 个大小为 \frac{a+b}{\gcd(a,b)} 的环,记 g=\gcd(a,b),n=\frac{a+b}{\gcd(a,b)},根据 Burnside 引理,等价类数量 = 所有置换下不动点数量的平均值,即

|X/G|=\frac{1}{|G|}\sum_{g\in G}X^g

其中 X^g 是在置换 g 的作用下的不动点数量。

在本问题中,|G|=nG 中的所有置换是「把 g 个大小为 n 的环旋转 i 个单位」,不动点数量是 (2^{\gcd(n,i)})^g=2^{g\gcd(n,i)}

总的答案是

ans=\sum_{i=1}^n2^{g\gcd(n,i)}

开始莫比乌斯反演

\begin{aligned} ans&=\sum_{d|n}2^{gd}\sum_{i=1}^\frac{n}{d}[\gcd(i,\frac{n}{d})=1]\\ &=\sum_{d|n}2^{gd}\varphi(\frac{n}{d}) \end{aligned}

预处理 \varphi 以及所有数的约数即可,复杂度 O((a+b)\ln(a+b)+T\cdot B),其中 B1\sim a+b 中约数个数最多的数的约数个数,为 240(在 a+b=720720 时取得)。

#include <bits/stdc++.h>
using namespace std;
const int N=1e6+5,mod=1000003;
inline int quick_pow(int a,int b){
    int ret=1;
    for(;b;b>>=1,a=1ll*a*a%mod)if(b&1)ret=1ll*ret*a%mod;
    return ret;
}
inline int dec(int a,int b){return (a-b<0)?a-b+mod:a-b;}
inline int Inv(int a){return quick_pow(a,mod-2);}
int phi[N],head[N],prime[N],pcnt,Pow[N];
bool vis[N];
struct Edge{
    int val,next;
}fac[N*14];
int tot=0;
inline void adde(int u,int v){
    fac[++tot]=(Edge){v,head[u]};head[u]=tot;
}
inline void init(int n){
    Pow[0]=1;
    for(int i=1;i<=n;++i){
        Pow[i]=2ll*Pow[i-1]%mod;
        for(int j=i;j<=n;j+=i)
            adde(j,i);
    }
    phi[1]=1;
    for(int i=2;i<=n;++i){
        if(!vis[i]){
            prime[++pcnt]=i;
            phi[i]=i-1;
        }
        for(int j=1;j<=pcnt&&prime[j]*i<=n;++j){
            vis[i*prime[j]]=1;
            if(i%prime[j]==0){
                phi[i*prime[j]]=phi[i]*prime[j];
                break;
            }
            phi[i*prime[j]]=phi[i]*(prime[j]-1);
        }
    }
}
int a,b,T;
int main(){
    init(1000000);
    scanf("%d",&T);
    while(T--){
        scanf("%d %d",&a,&b);
        if(!a||!b){
            puts("0");
            continue;
        }
        int g=__gcd(a,b),n=(a+b)/g,ret=0;
        for(int i=head[n],d;i;i=fac[i].next){
            d=fac[i].val;
            ret=(ret+1ll*Pow[d*g]*phi[n/d])%mod;
        }
        printf("%lld\n",dec(Pow[a+b],1ll*ret*g%mod*Inv(a+b)%mod));
    }
    return 0;
}