BZOJ 3601

Description:

给定一个非负整数 和正整数 ,令 为所有小于 且与 互质的数的 次方和。 对于给定的 的值。
由于 很大,所以给出 的质因数分解式。
$$N=\prod_{i=1}^{w} p_i^{a_i}(2<=p_i<=10^9~,~1<=a_i<=10^9)$$



Solution:

莫比乌斯反演之后有:

显然,可以构造数组 使得

那么假设我们构造出了数组 ,那么

显然这是个积性函数, 所以

所以 ,这样可以轻松做到 的时间复杂度。但是 怎么求呢,这很简单,高斯消元即可。所以总时间复杂度


Code:

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <iostream>
 
using namespace std;
 
const long long Mod=1e9+7;
 
long long fc[110][110]={{0}};
long long a[110]={0};
int d,W;
long long p_N[110]={0};
int D[1010]={0};
long long temp[1010]={0};
 
long long power(long long a,long long k)
{
    long long o=1;
    if(k<0) k+=Mod-1;
    for(;k>0;k>>=1)
    {
        if(k&1) o=o*a%Mod;
        a=a*a%Mod;
    }
    return o;
}
 
void Gauss()
{
    for(int i=1;i<=d+2;i++)
    {
        int g=i;
        for(int j=i;j<=d+2;j++)
            if(fc[j][i]!=0)
            {
                g=j;
                break;
            }
        if(g!=i) swap(fc[i],fc[g]);
        for(int j=1;j<=d+2;j++)
            if(j!=i && fc[j][i]!=0)
            {
                long long tmp=fc[j][i]*power(fc[i][i],Mod-2)%Mod;
                for(int p=0;p<=d+2;p++)
                    fc[j][p]=(fc[j][p]-tmp*fc[i][p]%Mod+Mod)%Mod;
            }
    }
    for(int i=1;i<=d+2;i++)
        a[i-1]=fc[i][0]*power(fc[i][i],Mod-2)%Mod;
    return;
}
 
long long N=1;
 
int main()
{
    scanf("%d%d",&d,&W);
    long long sum=0;
    for(int i=1;i<=d+2;i++)
    {
        sum=(sum+power(i,d))%Mod;
        fc[i][0]=sum;
        for(int j=1;j<=d+2;j++)
            fc[i][j]=power(i,j-1);
    }
    Gauss();
    for(int i=1;i<=W;i++)
    {
        int p,q;
        scanf("%d%d",&p,&q);
        N=N*power(p,q)%Mod;
        D[i]=p;
        temp[i]=power(p,Mod-2);
    }
    p_N[0]=1;
    for(int i=1;i<=d+1;i++)
        p_N[i]=p_N[i-1]*N%Mod;
    long long ans=0;
    for(int i=d+1;i>=0;i--)
    {
        long long cnt=1;
        for(int j=1;j<=W;j++)
        {
            cnt=cnt*(1-temp[j]+Mod)%Mod;
            temp[j]=temp[j]*D[j]%Mod;
        }
        ans=(ans+a[i]*p_N[i]%Mod*cnt)%Mod;
    }
    cout<<ans<<endl;
    return 0;
}