loj #2000. 「SDOI2017」数字表格 （莫比乌斯）

题目链接：https://loj.ac/problem/2000

代码：

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int MAXN=1e6+5;
const int MOD=1e9+7;
bool check[MAXN+10]; 
int prime[MAXN+10],f[MAXN][2],sum[MAXN],in[MAXN];
int mu[MAXN+10];
ll qpow(ll a,ll b)
{
	ll ans=1;a%=MOD;
    for(ll i=b;i;i>>=1,a=a*a%MOD)
        if(i&1)ans=ans*a%MOD;
    return ans;
}
ll inv(ll x)
{
	return qpow(x,MOD-2);
}
void Moblus()
{
	memset(check,false,sizeof(check));
	mu[1]=1;sum[1]=1;
	int tot=0;
	f[1][0]=1,f[0][0]=0;
	f[1][1]=1,f[0][1]=1;
	for(int i=2;i<MAXN;i++)
	{
		f[i][0]=(f[i-1][0]+f[i-2][0])%MOD;
		f[i][1]=inv(f[i][0]);
	}
	for(int i=2;i<=MAXN;i++)
	{
		sum[i]=1;
		if(!check[i])
		{
			prime[tot++]=i;
			mu[i]=-1;
		}
		for(int j=0;j<tot;j++)
		{
			if(i*prime[j]>MAXN) 
				break;
			check[i*prime[j]]=true;
			if(i%prime[j]==0)
			{
				mu[i*prime[j]]=0;
				break;
			}
			else
			{
				mu[i*prime[j]]=-mu[i];
			}
		}
	}
	sum[0]=in[0]=1;
	for(int i=1;i<MAXN;i++)
	{
		for(int j=i;j<MAXN;j+=i)
		{
			if(mu[j/i]==-1)
				sum[j]=1LL*sum[j]*f[i][1]%MOD;
			if(mu[j/i]==1)
				sum[j]=1LL*sum[j]*f[i][0]%MOD;
		}
	}
	for(int i=1;i<MAXN;i++)
	{
		sum[i]=1LL*sum[i-1]*sum[i]%MOD;
		in[i]=inv(sum[i]);
	}
}
ll cal(int n,int m)
{
	ll ret=1;
	int last;
	if(n<m)
		swap(n,m);
	for(int i=1;i<=m;i=last+1)
	{
		last=min(n/(n/i),m/(m/i));
		ret=1LL*ret*qpow(1LL*sum[last]*in[i-1]%MOD,1LL*(n/i)*(m/i))%MOD;
	}
	return ret;
}
int main()
{
	//freopen("in.txt","r",stdin);
	//freopen("out.txt","w",stdout);
	Moblus();
	int T;
	scanf("%d",&T);
	while(T--)
	{
		int n,m;
		scanf("%d%d",&n,&m);
		printf("%lld\n",cal(n,m));
	}
	return 0;
}