UVa 106/POJ 1305 Fermat vs. Pythagoras(数论&勾股数)

106 - Fermat vs. Pythagoras

Time limit: 3.000 seconds

http://uva.onlinejudge.org/index.php?option=com_onlinejudge&Itemid=8&category=115&page=show_problem&problem=42

http://poj.org/problem?id=1305

题意：

求≤N的本原勾股数组(PPT)的组数和≤N的非勾股数的个数。

思路：

参考《数论概论》P10的如下定理：

则对所有满足3≤s≤√(2n-1)，1≤t≤max{ s-2 , √(2n-s^2) }，gcd(s,t)=1的奇数s，t，求PPT和对应范围内的派生勾股数（本原勾股数组(a,b,c)的派生勾股数组为(ka,kb,kc)，k>1），然后标记已出现过的数。

最后统计未标记的数的个数即可。

复杂度：O(NlogN)

完整代码：

/*UVa: 0.192s*/
/*POJ: 0ms,1136KB*/

#include <cstdio>
#include <cstring>
#include <cmath>
const int maxn = 1000001;

bool vis[maxn];

int gcd(int a, int b) {return b ? gcd(b, a % b) : a;}

int main()
{
	int n, s, t, maxs, maxt, i, a, b, c;
	while (~scanf("%d", &n))
	{
		int cnt = 0;//PPT的组数
		memset(vis, 0, sizeof(vis));

		//循环生成所有的s和t，根据公式，s的最大值maxs<=√(2n-1)
		maxs = (int)sqrt((double)((n << 1) - 1));
		for (s = 3; s <= maxs; s += 2)
		{
			//求出t的最大值maxt，maxt<=s-2 && maxt<=√(2n-s^2)
			maxt = (int)sqrt((double)((n << 1) - s * s));
			if (s - 2 < maxt) maxt = s - 2;
			for (t = 1; t <= maxt; t += 2)
			{
				if (gcd(s, t) == 1)
				{
					++cnt;
					a = s * t;
					b = (s * s - t * t) >> 1;
					c = (s * s + t * t) >> 1;

					//标记PPT及其派生勾股数
					for (i = 1; c * i <= n; ++i)
						vis[a * i] = vis[b * i] = vis[c * i] = true;
				}
			}
		}
		printf("%d ", cnt);

		//统计没有出现过的数的个数
		cnt = 0;
		for (i = 1; i <= n; ++i)
			if (!vis[i]) ++cnt;
		printf("%d\n", cnt);
	}
	return 0;
}