SDOI2017 数字表格

本文移植自个人的原博客（原文）

题目描述

题解

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
using namespace std;
typedef long long ull;
typedef double db;;
const int mod = 1000000007;
const int N = 1e6+1;
ull fi[N],p[N],miu[N],phi[N],tot,d[N],f[N],fp[N][3],ff[N],fs[N],rev[N];
bool vis[N];
int ca1[33],ca2[33];

ull qpow(ull a,ull b)
{
    ull ret=1;
    if(b<0)b+=mod-1;
    while(b)
    {
        if(b&1)ret=(ret*a)%mod;
        b>>=1;
        a=(a*a)%mod;
    }
    return ret;
}
void init()
{
    miu[1]=phi[1]=fs[0]=f[1]=1;
    for(int i=2;i<N;++i)
    {
        f[i]=(f[i-1]+f[i-2])%mod;
        if(!vis[i])miu[i]=-1,p[++tot]=i,phi[i]=i-1;
        for(int j=1;j<=tot&&p[j]*i<N;++j)
        {
            vis[i*p[j]]=true;
            if(i%p[j]==0){miu[i*p[j]]=0;phi[i*p[j]]=phi[i]*p[j];break;}
            miu[i*p[j]]=-miu[i];
            phi[i*p[j]]=phi[i]*(p[j]-1);
        }
    }

    for(int i=1;i<N;++i)
    fp[i][0]=qpow(f[i],-1),fp[i][1]=1,fp[i][2]=f[i];

    for(int i=1;i<N;++i)ff[i]=1;

    for(int i=1;i<N;++i)
    for(int j=i;j<N;j+=i)
    ff[j]=(ull)ff[j]*fp[i][miu[j/i]+1]%mod;
    for(int i=1;i<N;++i)fs[i]=(ull)fs[i-1]*ff[i]%mod;

    rev[N-1]=qpow(fs[N-1],-1);

    for(int i=N-2;i>=0;--i)rev[i]=(ull)rev[i+1]*ff[i+1]%mod;
    return;
}
//ull solve(int n,int m)
//{
//    ull ret=0;
//    int top=min(n,m);
//    ull i=1,j;
//    while(i<=top)
//    {
//        j=min(n/(n/i),m/(m/i));
//        (ret+=(ull)(miu[j]-miu[i-1])*((n/i)*(m/i))%mod)%=mod; 
//        i=j+1;
//    }
//    return ret;
//}
//ull mobi(ull n,ull m)
//{
//    if(n>m)swap(n,m);
//    ull ans=0;
//    for (int i=1,r;i<=n;i=r+1)
//    {
//        r=min(n/(n/i),m/(m/i));
//        ans=(ans+sum(n/i)*sum(m/i)%mod*(miu[r]-miu[i-1])%mod)%mod;
//    }
//    return ans;
//}
//ull cheng(ull a,ull b)
//{
//    ull ret=0;
//    int tmp1=0,tmp2=0;
//    memset(ca1,0,sizeof(ca1));
//    memset(ca2,0,sizeof(ca2));
//    while(a)ca1[++tmp1]=a%10,a/=10;
//    while(b)ca2[++tmp2]=b%10,b/=10;
//    int *x=ca1,*y=ca2;
//    if(tmp2>tmp1)swap(x,y),swap(tmp1,tmp2);
//    for(int i=1;i<=tmp2;++i)x[i]*=y[i];
//    for(int i=1;i<=tmp2;++i)
//    if(x[i]>10)
//    {
//        x[i+1]+=x[i]/10;
//        x[i]%=10;
//    }
//    for(int i=1;i<=tmp1;++i)y[tmp1+1-i]=x[i];
//    for(int i=1;i<=tmp1;++i)ret=(ret*10+(y[i]))%mod;
//    return ret;
//}
int main()
{
    int T,n,m,k;
//    freopen("c.out","w",stdout);
    init();
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d",&n,&m);
        ull ans=1;
        int mn=min(n,m);
        for(int i=1,r;i<=mn;i=r+1)
        {
            r=min(n/(n/i),m/(m/i));
            ans=ans*qpow(fs[r]*rev[i-1]%mod,(ull)(n/i)*(m/i)%(mod-1))%mod;
        }
        cout<<ans<<endl;
    }
}