loj SDOI2017数字表格

太才了

\(\prod \limits_{i=1}^{n}\prod\limits_{j=1}^{m}f[gcd(i,j)]\)
\(\prod\limits_{k=1}^{n}f[k]^{\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[gcd(i,j)==k]}\)
幂我们很熟悉
就是
\(g(x)=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[gcd(i,j)==x]=\sum\limits_{i=1}^{\frac{n}{x}}\sum\limits_{j=1}^{\frac{m}{x}}[gcd(i,j)==1]\)
\(f(x)=\sum\limits_{x|d}^ng(d)=\frac{n}{x}\frac{m}{x}\)
\(g(x)=\sum\limits_{x|d}^n\mu(\frac{d}{x})f(d)\)
\(g(x)=\sum\limits_{x|d}^n\mu(\frac{d}{x})\frac{n}{d}\frac{m}{d}\)
\(g(x)=\sum\limits_{d=1}^{n/x}\mu(d)\frac{n}{xd}\frac{m}{xd}\)
带回去
\(\prod\limits_{k=1}^{n}f[k]^{\sum\limits_{d=1}^{n/k}\mu(d)\frac{n}{dk}\frac{m}{dk}}\)
令i=d*k
\(\prod\limits_{k=1}^{n}f[k]^{\sum\limits_{k|i}^{n}\mu(i/k)\frac{n}{i}\frac{m}{i}}\)
\(\prod\limits_{k=1}^{n}\prod\limits_{k|i}^{n}f[k]^{\mu(i/k)\frac{n}{i}\frac{m}{i}}\)
\(\prod\limits_{i=1}^{n}\prod\limits_{k|i}^{n}f[k]^{\mu(i/k)\frac{n}{i}\frac{m}{i}}\)
\(Au(k)=\prod\limits_{k|i}^{n}f[k]^{\mu(i/k)}\)
\(ans=\prod\limits_{i=1}^{n}Au(k)^{\frac{n}{i}\frac{m}{i}}\)

注意,错误

幂次是不能直接模除的,用\(x^{p-1}\equiv 1(mod p)\)

代码

#include <bits/stdc++.h>
#define ll long long
using namespace std;
const ll N=1e6+7,mod=1e9+7;
ll read() {
    ll x=0,f=1;char s=getchar();
    for(;s>'9'||s<'0';s=getchar()) if(s=='-') f=-1;
    for(;s>='0'&&s<='9';s=getchar()) x=x*10+s-'0';
    return x*f;
}
ll f[N],mu[N],pri[N],Au[N],cnt;
bool vis[N];
ll q_pow(ll a,ll b) {
	ll ans=1;
	while(b) {
		if(b&1) ans=ans*a%mod;
		a=a*a%mod;
		b>>=1;
	}
	return ans;
}
void Euler() {
	mu[1]=vis[1]=1;
	for(ll i=2;i<=1000000;++i) {
		if(!vis[i]) {pri[++cnt]=i,mu[i]=-1;}
		for(ll j=1;j<=cnt&&i*pri[j]<=1000000;++j) {
			vis[i*pri[j]]=1;
			if(i%pri[j]==0) {
				mu[i*pri[j]]=0;
				break;
			} else mu[i*pri[j]]=-mu[i];
		}
	}
	f[1]=1;
	for(ll i=2;i<=1000000;++i) f[i]=(f[i-1]+f[i-2])%mod;
	for(ll i=0;i<=1000000;++i) Au[i]=1;
	for(ll i=1;i<=1000000;++i) {
		ll inv=q_pow(f[i],mod-2);
		for(ll j=i;j<=1000000;j+=i) {
			if(mu[j/i]==-1) Au[j]=Au[j]*inv%mod;
			else if(mu[j/i]==1) Au[j]=Au[j]*f[i]%mod;
		}
	}
	for(ll i=1;i<=1000000;++i) Au[i]=(Au[i]*Au[i-1])%mod;
}
int main() {
	Euler();
	ll T=read();
	while(T--) {
		ll n=read(),m=read();
		if(n>m) swap(n,m);
		ll ans=1;
		for(ll l=1,r=1;l<=n;l=r+1) {
			r=min(n/(n/l),m/(m/l));
			ans=ans*q_pow(Au[r]*q_pow(Au[l-1],mod-2)%mod,(n/l)*(m/l)%(mod-1))%mod;
		}
		cout<<ans<<"\n";
	}
    return 0;
}

posted @ 2019-02-15 12:35  ComplexPug  阅读(322)  评论(0编辑  收藏  举报