POJ 2480 Longge's problem

#include<cstdio>
using namespace std;
typedef long long ll;
ll p[6010],prime[50010];
ll l,d,ans,s[20],k[20],b[20];
void get_prime(){
	ll i,j;
	prime[0]=prime[1]=1;
	for(i=2;i<=50000;i++)
		if(!prime[i]){
			if(i<245)
				for(j=i*i;j<=50000;j+=i)
					prime[j]=1;
			p[++l]=i;
		}
}
ll jc(ll n,ll m){
	ll i,sum=1;
	for(i=1;i<=m;i++)
		sum*=n;
	return sum;
}
ll euler(ll n,ll m){
	if(n==1||m==0)return 1;
	return jc(n,m)/n*(n-1);
}
void search(ll n,ll m,ll x){
	ll i,sum;
	if(n>d){
		sum=(x/m);
		for(i=1;i<n;i++)
			sum*=euler(s[i],b[i]);
		ans+=sum;
		return;
	}
	for(i=0;i<=k[n];i++){
		b[n]=i;
		search(n+1,m*jc(s[n],i),x);
	}
}
ll fj(ll n){
	ll i,x=n;
	d=0;
	for(i=1;i<=l&&n!=1;i++)
		if(n%p[i]==0){
			s[++d]=p[i];
			k[d]=0;
			while(n%p[i]==0){
				n/=p[i];
				k[d]++;
			}
		}
	if(n!=1){
		s[++d]=n;
		k[d]=1;
	}
	ans=0;
	search(1,1,x);
	return ans;
}
int main(){
	ll n;
	get_prime();
	while(scanf("%lld",&n)!=EOF)
		printf("%lld\n",fj(n));
	return 0;
}