本来挺简单的,可是n太大了,用欧拉线性筛都不行啊!这可咋办呢?
P为N的质因数
Phi(N)=N*(1-(1-1/p1)) * (1-(1-1/p2))… * (1-(1-1/pn)).
然后呢?这么大的数让我分解个毛线!
分解大的数需要用Pollard-rho整数分解
实现方法:生成两个整数a和b,计算p=gcd(a-b,n),直到p不为1或者a,b出现循环为止,若p=n,则p为质数,否则p为n的一个约数。选取一个小的随机数x1,迭代生成xi=xi-1 ^2+k,一般取k=1,若序列出现循环则退出。计算p=gcd(xi-1 - xi,n),若p=1,返回上一步,直到p>1为止。若p=n,则n为素数,否则p为n的一个约数并递归分解p和n/p.
我的理解:随机生成n的因数
很大的数怎么判断是不是质数啊!
用Miller-Rabin!
费马小定理:对于素数p和任意整数a,有ap ≡ a(mod p)(同余)。反过来,满足ap ≡ a(mod p),p也几乎一定是素数。
伪素数:如果n是一个正整数,如果存在和n互素的正整数a满足 an-1 ≡ 1(mod n),我们说n是基于a的伪素数。如果一个数是伪素数,那么它几乎肯定是素数。
Miller-Rabin测试:不断选取不超过n-1的基b(s次),计算是否每次都有bn-1 ≡ 1(mod n),若每次都成立则n是素数,否则为合数.
注意,MIller-Rabin测试是概率型的,不是确定型的,不过由于多次运行后出错的概率非常小,所以实际应用还是可行的。(一次Miller-Rabin测试其成功的概率为3/4)
下面还有一个定理,能提高Miller测试的效率:
二次探测定理:
如果p是奇素数,则 x2 ≡ 1(mod p)的解为 x = 1 || x = p - 1(mod p);
注意啦:我们在进行素数测试时,mod可能很大,所以要化乘为加(类似于快速幂算法)。
#include <cstdio>
#include <iostream>
#include <cstdlib>
#include <algorithm>
#define ll long long
using namespace std;
ll prime[99999],cnt;
ll gcd(ll a,ll b)
{
if(!b) return a;
return gcd(b,a%b);
}
ll mul(ll x,ll y,ll mod)//快速加
{
ll ans=0;
while(y)
{
if(y%(1ll*2)) ans=(ans+x%mod)%mod;
y/=(1ll*2);
x=(x%mod+x%mod)%mod;
}
return ans;
}
ll fast_pow(ll x,ll y,ll m)//快速幂
{
ll ans=1;
x%=m;
while(y)
{
if(y%2)
ans=mul(ans, x, m);
y/=2;
x=mul(x, x, m);
}
return ans;
}
bool MR(ll n)
{
if(n == 2) return 1;
if(n < 2 || !(n & 1)) return 0;
ll m=n-1;
ll k=0;
while((k&1)==0)
{
k++;
m/=2;
}
for(int i=0; i<12; i++)
{
ll a=rand()%(n - 1)+1;
ll x=fast_pow(a, m, n);
ll y=0;
for(int j=0; j<k; j++)
{
y=mul(x,x,n);
if(y==1&&x!=1&&x!=n-1) return 0;
x = y;
}
if(y!=1) return 0;
}
return 1;
}
ll rho(ll n,ll c)//找因子
{
ll i=1,k=2;
ll x=rand()%(n-1)+1;
ll y=x;
while(1)
{
i++;
x=(mul(x,x,n)+c)%n;
ll d=gcd(((y-x)+n)%n,n)%n;
if(d>1&&d<n) return d;
if(y==x) return n;
if(i==k)//一个优化,我TM也不知道为啥
{
y=x;
k<<=1;
}
}
}
void find(ll n,ll c)//分解的递归过程
{
if(n==1) return;
if(MR(n))
{
prime[++cnt]=n;
return;
}
ll p=n;
ll k=c;
while(p>=n) p=rho(n,c--);
find(n/p,c);
find(p,c);
}
int main()
{
ll n;
scanf("%lld",&n);
find(n,120);
sort(prime+1,prime+cnt+1);
ll t=unique(prime+1,prime+cnt+1)-prime-1;//去重函数,不要忘了去了重复的质数啊
/* for(int i=1;i<=t;i++)
printf("%lld ",prime[i]);*/
double ans=n;
for(int i=1;i<=t;i++)
ans*=(1.0-(1.0/(double)prime[i]));
printf("%lld",(ll) ans);
}