[2012集训队互测]JZPKIL - 生成函数,伯努利数,数论,莫比乌斯反演,狄利克雷卷积

Pollard_Rho质因数分解  似乎有点卡常？

#include"bits/stdc++.h"

using namespace std;

const int N=3005,P=1000000007;
typedef long long LL;

int inv(LL a,LL t=P-2){
	int r=1;t%=P-1,a%=P;
	if(t<0)t+=P-1;
	while(t){
		if(t&1)r=1LL*r*a%P;
		a=1LL*a*a%P;t>>=1;
	}
	return r;
}

int b[N],f[N][N],val[N][N],C[N][N];
void pre(){
	b[0]=C[0][0]=1;
	for(int i=1;i<=3001;i++)for(int j=0;j<=i;j++)C[i][j]=(C[i-1][j-1]+C[i-1][j])%P;
	for(int d=1;d<=3000;d++){
		for(int k=0;k<d;k++)b[d]=(b[d]-1LL*b[k]*C[d+1][k]%P+P)%P;
		b[d]=1LL*b[d]*inv(d+1)%P; 
	} 
	for(int d=0;d<=3000;d++){
		int inv_d=inv(d+1);
		for(int k=0;k<=d;k++)val[d][d+1-k]=1LL*b[k]*C[d+1][k]%P*inv_d%P;
		(val[d][d]+=1)%=P;
	}
}

LL mult(LL x,LL y,LL z){
	//__int128 r=x;r=r*y;
	//return (LL)(r%z);
	return (x*y-(LL)(((long double)x*y+0.5)/(long double)z)*z+z)%z;
}
//LL mult(LL a,LL b,LL c){return (LL)((long double)a*b-(LL)(((long double)a*b+0.5)/c)*(long double)c);}

LL ksm(LL a,LL t,LL P){
	LL r=1;
	while(t){
		if(t&1)r=mult(r,a,P);
		a=mult(a,a,P);t>>=1;
	}
	return r;
}

const int JudgeP[] = {2,3,5,7,11,13,17,19,23};

inline int irand(){return (rand()<<16)^rand();}

bool Miller_Rabin(LL n){ 
	if(n==2||n==3||n==5||n==7||n==11||n==13||n==17||n==19||n==23)return true;
	if((n&1)==0)return false;
	LL r=n-1;int k=0;
	while((r&1)==0)r>>=1,k++;
	for(int i=0;i<9;i++){
		LL x=ksm(JudgeP[i],r,n),y;
		for(int i=0;i<k;i++){
			y=mult(x,x,n);
			if(y==1&&x!=1&&x!=n-1)return false;
			x=y;
		}
		if(y!=1)return false;
	}
	return true;
}

int T,numc;
LL num[70];

LL gcd(LL x,LL y){return y?gcd(y,x%y):x;}

#define Func(x) (mult(x,x,n)+1) 

LL getnum(LL n){
	LL x=irand()%n+1,y=Func(x);
	while(x!=y){
		LL k=gcd(n,(y-x+n)%n);
		if(k!=1&&k!=n)return k;
		x=Func(x),y=Func(Func(y));
	}
	return -1;
}

void Pollard_Rho(LL n){
	while(n%2==0){
		num[++numc]=2;
		n>>=1;
	}
	if(n==1)return;
	if(Miller_Rabin(n)){
		num[++numc]=n;
		return;
	}
	LL x;
	do x=getnum(n);while(x==-1);
	Pollard_Rho(x);
	Pollard_Rho(n/x);
}

LL n;

void work(){
	srand(19260817);
	scanf("%d",&T);
	int x,y;
	while(T--){
		scanf("%lld%d%d",&n,&x,&y);
		numc=0;
		if(n==1){puts("1");continue;}
		Pollard_Rho(n);
		sort(num+1,num+numc+1);
		int l=unique(num+1,num+numc+1)-num-1;
		LL r=0,w,m,s,k,g;
		for(int i=1;i<=y+1;i++){
			w=1,m=n;
			for(int j=1;j<=l;j++){
				s=1,k=0;
				while(m%num[j]==0)m/=num[j],s*=num[j],k++;
				g=0;
				for(int t=0;t<=k;t++)g=(g+inv(num[j],1LL*t*(x-i)+1LL*k*i))%P;
				for(int t=0;t<k;t++)g=(g-inv(num[j],1LL*t*(x-i)+1LL*(k-1)*i+y))%P;
				w=1LL*w*(g+P)%P;
			}
			r=((r+1LL*val[y][i]*w)%P+P)%P;
		}
		r=1LL*r*inv(n,y)%P;
		printf("%lld\n",r);
	}
}

int main(){
	pre();
	work();
	return 0;
}