累乘函数线性逆元打表，阶乘反演——bzoj4816

学了一种新套路，倒序打表函数的逆元可以直接线性完成

#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define LL long long
const long long mod = 1e9+7;
#define maxn 1000001

LL Pow(LL a,LL b){
    if (!b)
        return 1LL;
    LL x=Pow(a,b/2);
    x=x*x%mod;
    if (b&1LL)
        x=x*a%mod;
    return x;
}
/*ll Pow(ll a,ll b){
    ll res=1;
    while(b){
        if(b%2)
            res=res*a%mod;
        b>>=1;a=a*a%mod;
    }
    return res;
}*/
ll n,m;
ll F[maxn],pre[maxn],invF[maxn];
void init1(){
    F[0]=0;F[1]=1;
    for(int i=2;i<maxn;i++)
        F[i]=(F[i-1]+F[i-2])%mod;
    pre[0]=1;
    for(int i=1;i<maxn;i++)//临时数组算累乘 
        pre[i]=pre[i-1]*F[i]%mod;

    ll tmp=Pow(pre[maxn-1],mod-2);
    for(int i=maxn-1;i>=1;i--)
        invF[i]=tmp*pre[i-1]%mod,tmp=tmp*F[i]%mod;
}

bool vis[maxn];
ll prime[maxn],mm,mu[maxn];
void init2(){
    mu[1]=1; 
    for(int i=2;i<maxn;i++){
        if(!vis[i]){
            mu[i]=-1;
            prime[++mm]=i;
        }
        for(int j=1;j<=mm;j++){
            if(i*prime[j]>=maxn)break;
            vis[i*prime[j]]=1;
            if(i%prime[j]==0){
                mu[i*prime[j]]=0;
                break;
            }
            else mu[i*prime[j]]=-mu[i];
        }
    }
}

ll mul[maxn],invmul[maxn],g[maxn];
void init3(){
    for(int i=1;i<maxn;i++)g[i]=1;
    for(int i=1;i<maxn;i++)
        for(int j=1;j*i<maxn;j++){
            if(mu[j]==1)
                g[i*j]=g[i*j]*F[i]%mod;
            else if(mu[j]==-1)
                g[i*j]=g[i*j]*invF[i]%mod;
        }
    mul[0]=1;
    for(int i=1;i<maxn;i++)    
        mul[i]=mul[i-1]*g[i]%mod;
    invmul[maxn-1]=Pow(mul[maxn-1],mod-2);
    for(int i=maxn-2;i>=0;i--)
        invmul[i]=invmul[i+1]*g[i+1]%mod;
}

int main(){
    int t;cin>>t;
    init1();
    init2();
    init3();
    while(t--){
        cin>>n>>m;
        if(n>m)swap(n,m);
        ll ans=1;
        for(int l=1,r;l<=n;l=r+1){
            r=min(n/(n/l),m/(m/l));
            ll tmp=mul[r]*invmul[l-1]%mod;
            ans=ans*Pow(tmp,(n/l)*(m/l)%(mod-1))%mod;
        }
        cout<<ans<<'
';
    }
}