JZOJ5787 轨道 的代码

JZOJ5787 轨道 的代码

#include<bits/stdc++.h>
using namespace std;
const int maxn=30001;
const int maxm=3200;
const int maxv=1e7;
const int mod=10007;
int n,m,K;
int dp[maxn+5][maxm+5];
int a[maxm+5],at;
int mp[maxv+5];
int mu[maxv+5];
bool vis[maxv+5];
int pr[maxv+5],prt;
int V[maxm+5][maxm+5];
signed main(){
	int t;
	scanf("%d%d%d",&n,&m,&K);
	mu[1]=1;
	for(int i=2;i<=K;i++){
		if(!vis[i]) pr[++prt]=i,mu[i]=-1; 
		for(int j=1;j<=prt&&pr[j]*i<=K;j++){
			vis[pr[j]*i]=1;
			if(!(i%pr[j])){mu[pr[j]*i]=0;break;}
			mu[pr[j]*i]=-mu[i];
		}
	}
	for(int i=1;i*i<=K;i++){
		if(K%i==0){
			a[++at]=i;
			if(i*i!=K) a[++at]=K/i;
		}
	}
	sort(a+1,a+at+1);
	for(int i=1;i<=at;i++){
		mp[a[i]]=i;
		t=m/a[i];
		for(int j=1;j<=at;j++)
			dp[1][i]=(dp[1][i]+mu[a[j]]*(t/a[j]))%mod;
	}
	for(int i=1;i<=at;i++)
		for(int j=1;j<=i;j++)
			if(a[i]%a[j]==0)
				V[i][++V[i][0]]=j;
	for(int i=2;i<=n;i++){
		for(int j=1;j<=at;j++){
			for(int k=1;k<=V[j][0];k++){
				t=V[j][k];
				dp[i][j]=(dp[i][j]+dp[i-1][t]
					*dp[1][mp[a[j]/a[t]]])%mod;
			}
		}
	}
	printf("%d",(dp[n][at]+mod)%mod);
	return 0;
}