「SDOI2017」数字表格

「SDOI2017」数字表格

链接:https://www.luogu.com.cn/problem/P3704


显然答案是:

i=1nj=1mfgcd(i,j)

化简考虑枚举 f 的下标:

=di=1ndj=1mdfd[gcd(i,j)=1]=dfdi=1ndj=1mdk|gcd(i,j)μ(k)=dfdk=1min(ndmd)μ(k)×ndkmdk

这个式子我们可以外层数论分块内层数论分块在 O(n) 的时间内算出,但是显然不够。继续优化式子。

T=k×d。将 T 提出,枚举 T

那么原式为:

T=1(d|Tfdμ(Td))nTmT

而对于 d|Tfdμ(Td) 我们可以枚举倍数在 O(nlnn) 的时间内求出。

外层再用数论分块优化,算上求逆元,单次查询是 nlogMOD 的。

总复杂度是 O(TnlogMOD+nlnnlogMOD)

代码如下:

#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int MOD = 1e9+7;
const int MAXN = 1e6+5;
const int MX = 1e6;
int n,m,prime[MAXN],p[MAXN],pc;
ll f[MAXN],Sum[MAXN],mu[MAXN],g[MAXN],inv[MAXN];
ll qpw(ll x,ll b)
{
	ll now=1,tmp=x,ans=1;
	while(now<=b)
	{
		if(now&b) ans=ans*tmp%MOD;
		tmp=tmp*tmp%MOD;
		now<<=1;
	}
	return ans;
}
int main()
{
	int T;scanf("%d",&T);
	f[0]=0;f[1]=1;
	for(int i=2;i<=MX;++i) f[i]=(f[i-1]+f[i-2])%MOD;
	for(int i=0;i<=MX;++i) inv[i]=qpw(f[i],MOD-2);
	mu[1]=1;
	for(int i=2;i<=MX;++i)
	{
		if(!p[i]) prime[++pc]=i,mu[i]=-1;
		for(int j=1;j<=pc&&1ll*prime[j]*i<=MX;++j)
		{
			p[i*prime[j]]=1;
			if(i%prime[j]==0)
			{
				mu[i*prime[j]]=0;
				break;
			}
			mu[i*prime[j]]=-mu[i];
		}
	}
	for(int i=0;i<=MX;++i) g[i]=1;
	for(int i=1;i<=MX;++i)
		for(int j=i;j<=MX;j+=i)
			if(mu[j/i]!=0) g[j]=g[j]*(mu[j/i]==1?f[i]:inv[i])%MOD;
	for(int i=2;i<=MX;++i) g[i]=g[i]*g[i-1]%MOD;
	while(T--)
	{
		scanf("%d %d",&n,&m);
		ll Ans=1;
		for(ll l=1,r;l<=min(n,m);l=r+1)
		{
			r=min(n/(n/l),m/(m/l));
			Ans=Ans*qpw(g[r]*qpw(g[l-1],MOD-2)%MOD,(n/l)*(m/l))%MOD;
		}
		printf("%lld\n",Ans);
	}
	return 0;
}
posted @   夜空之星  阅读(45)  评论(0编辑  收藏  举报
相关博文:
阅读排行:
· 地球OL攻略 —— 某应届生求职总结
· 周边上新:园子的第一款马克杯温暖上架
· Open-Sora 2.0 重磅开源!
· 提示词工程——AI应用必不可少的技术
· .NET周刊【3月第1期 2025-03-02】
点击右上角即可分享
微信分享提示