解题报告之 SOJ1839 Relatives 欧拉函数算法详解

最新推荐文章于 2018-07-26 00:09:27 发布

StanleyClinton

最新推荐文章于 2018-07-26 00:09:27 发布

阅读量1k

点赞数 1

分类专栏： ACM_数论 ACM 文章标签： SOJ1839 Relatives 欧拉函数素数筛线性筛法数论

本文链接：https://blog.csdn.net/maxichu/article/details/45540297

版权

ACM 同时被 2 个专栏收录

118 篇文章 3 订阅

订阅专栏

ACM_数论

20 篇文章 2 订阅

订阅专栏

解题报告之 SOJ1839 Relatives

Description

Problem

Given n, a positive integer, how many positive integers less than n are relatively prime to n? Two integers a and b are relatively prime if there are no integers x > 1, y > 0, z > 0 such that a = xy and b = xz.

There are several test cases. For each test case, standard input contains a line with n <= 1,000,000,000. A line containing 0 follows the last case.

For each test case there should be single line of output answering the question posed above.

Sample Input

7
12
0

Output for Sample Input

6
4

题目大意：给你一个数n，不超过1e9，让你找出不超过它且与之互质的数的个数。

分析：一道裸的欧拉函数，不过作为存模板的博文，裸的还是比较通用一些。一个数n 的不超过它的与之互诉的数记为 φ(n)。且

φ(n) = n * ( 1 - 1/p1 ) * ( 1 - 1/p2 ) * …… * ( 1 - 1/pm ) **

（其中pk表示n的某个不重复的质因子）

然后我们一般操作的方法是先令ans=n，一边分解n一边更新ans。每检测到一个质因子，则ans=ans*（1-1/pk）=ans-ans/pk。

最后将n分解完之后ans即是答案。此处我找质因子的方法是先素数筛然后用直接用素数去试验，比较高效。如果不明白素数筛的建议先参考我另一篇博文《素数判定Miller_Rabin 算法详解》中关于素数筛法的论述。因为后面我还要讲线性筛法批量求欧拉函数。

#include<iostream>
#include<cmath>
#include<algorithm>
using namespace std;

const int MAXN = 1e5;

int isprime[MAXN];
int prime[MAXN];
int cnt;
void getP()
{
	cnt = 0;
	for(int i = 1; i < MAXN; i++)
		isprime[i] = 1;
	for(int i = 2; i < MAXN; i++)
	{
		if(!isprime[i])continue;
		prime[cnt++] = i;
		for(int j = 2 * i; j < MAXN; j += i)
		{
			isprime[j] = 0;
		}
	}
}

int euler( int n )
{
	int ans = n;

	for(int i = 0; prime[i] * prime[i] <= n; i++)
	{
		if(n%prime[i] == 0)
		{
			ans = ans - ans / prime[i]; //ans=ans*(1-1/pi)
			while(n%prime[i] == 0)
			{
				n /= prime[i];
			}
		}
	}

	if(n > 1)
		ans = ans - ans / n;
	return ans;
}

int main()
{
	getP();
	int n;
	while(cin >> n&&n)
	{
		cout << euler( n ) << endl;
	}
}

---------------------------------------------------------------------------------------------------------------------------------------------

好了，这种方法适用于求单个数的欧拉函数且数特别大的时候，但是如果遇到需要求批量的欧拉函数的情况这种方法就显的不那么好用了。于是我们借鉴素数筛法，在筛素数的同时一并把欧拉函数也求出来。

首先，要明确一些性质，下面三条性质我们会在欧拉函数线性筛中使用到。

1.如果p为素数，则φ(p)=p-1；

2.若a为n的质因子，且(n%a==0&&(n/a)%a==0)，即如果a是n的多重质因子，则φ(n)=φ(n/a)*a；

3.若a为n的质因子，但(n%a==0&&(n/a)%a! =0)，即如果a是n的单质因子，则φ(n)=φ(n/a)*(a-1)；

我们首先保持素数筛的流程不变，再向里面添加一些功能。当我们在素数筛时，不再仅仅对素数更新，我们对每个数字都试图更新后面的部分数。如果发现一个未访问过的数为素数 i 时，加入到prime数组，并且让phi [ i ] =i-1；不管它是不是素数，我们都更新后面的数，我们遍历已经找到的素数prime数组，i*prime[j]表示更新目标，先将访问标志设为1，表示之后不用再更新了。接着此时prime[j]是它的质因子，即 i =（i*prime[j]）/prime[j]，可以用性质2或者性质3。

那么到底用性质2还是性质3呢？ -> 取决于a是多重还是单重质因子。 -> 怎么知道它是单重还是多重呢？ -> 看 i 中还能不能分解出一个prime[j]。-> 能分解则说明prime[j]是多重，不能分解则为单重。

void getPhi()
{
	phi[1] = 1;
	cnt = 0;
	for(int i = 1; i < MAXN; i++)
		isprime[i] = 1;

	for(int i = 2; i < MAXN; i++)
	{
		if(isprime[i])
		{
			prime[cnt++] = i;
			phi[i] = i - 1;
		}
		for(int j = 0; j<cnt&&prime[j]*i<MAXN; j ++)
		{
			isprime[i*prime[j]] = 0;
			if(i%prime[j] == 0)
			{
				phi[i*prime[j]] = phi[i] * (prime[j]);
				break;
			}
			else
			{
				phi[i*prime[j]] = phi[i]*( prime[j] - 1 );
			}
		}
	}
}

整个代码中最难理解的就是那个break了。这里的break是为了保证每个phi都只被更新一次，即它采用的策略是将每一个数通过它的最小质因子与剩下的因子来更新，所以一旦一个数字检测到能整除某个质数那么它更新了以当前质数prime[j]为最小因子且剩下因子为i的数，但是如果它继续更新prime[j+1]，那么i至少为prime[j]*p，那么注定在 i*=prime[j+1]*p时还会被prime[j]更新一次，那么就重复了。因此，这里要break，目的是为了保持所有数都是通过最小质因子和其剩余因子来更新的。

比如要更新12=2*2*3，当扫描到i=4时，prime[j]=2，此时更新8，并且发现4%2==0，那么就break，因为12会由最小质因子2和6来更新。下面论述一下这种策略如何保证所有数都被更新到且都只更新一次：

1.所有质数p在素数筛时会直接制定一个值p-1，所以都更新到了，且之后的更新都是两个数的乘积，不可能更新到素数；

2.对于合数来说，任意合数都能拆成其最小质因子和剩余因子积的形式，且这剩余因子积的phi在之前一定都更新过了，故所有合数通过这种策略来更新都更新到了且都只更新了一次。证毕。

---------------------------------------------------------------------------------------------------------------------------------------------

学习了线性筛，再提示一个进阶的知识吧，如何求一个数约数的个数，根据公式有n=p1^a1 * p2^a2 * …… * pn^an，那么这个数的约数个数 = (a1+1)*(a2+1)*(a3+1)*……*(an+1)。为什么是ak+1呢，表示选这个0~ak个该质因子，则有ak+1种选法。且约数函数具有可乘性。采用最小质因子*剩余因子积的策略，不难归纳出：

1.对于i和prime[j]，如果之前已经出现过最小质因子即i%prime[j]==0，那么约数个数divnum[i*prime[j]]=divnum[i]/(minp[i]+1)*(minp[i]+2)，minp表示数i的最小质因子数，为什么是/(minp[i]+1)*(minp[i]+2)呢？因为此时质因子多了一个，那么公式计算式对应乘数就应该+1。同时minp[i*prime[j]]=min[i]+1。注意同样要break.

2.如果未出现最小质因子，那么i%prime[j]！=0，那么根据可乘性divnum[i*prime[j]] = divnum[i] * divnum[prime[j]]; ，同时因为之前还没出现过最小质因子，那么minp[i*prime[j]]=1。

上代码：

void getDiv()
{
	for(int i = 1; i < MAXN; i++)
		isprime[i] = 1;
	for(int i = 2; i < MAXN; i++)
	{
		if(isprime[i])
		{
			prime[cnt++] = i;
			minp[i] = 1; //质数的最小质因子为其本身，数量为1
			divnum[i]= 2;//质数的约数为1和它本身
		}
		for(int j = 0; j < cnt&&i*prime[j] < MAXN; j++) //线性筛
		{
			isprime[i*prime[j]] = 0;
			if(i%prime[j] == 0)  //同样按照最小质因子*剩余因子积的策略
			{
				divnum[i*prime[j]] = divnum[i] / (minp[i] + 1)*(minp[i] + 2);	//根据公式进行更新
				minp[i*prime[j]] = minp[i] + 1; //最小质因子数多一个
				break;
			}
			else //如果之前没有最小质因子，说明prime[j]为最小质因子，且数量为1
			{
				divnum[i*prime[j]] = divnum[i] * divnum[prime[j]];	//约数函数的可乘性
				minp[i*prime[j]] = 1;
			}
		}
	}
}

嗯，这个应该比较全面了。