Pollard_rho大数质因数分解+拉格朗日四平方和定理（bzoj 2904: 平方和）

2904: 平方和

Time Limit: 10 Sec   Memory Limit: 128 MB
Submit: 160   Solved: 73
[ Submit][ Status][ Discuss]

Description

给定一个整数N，求N最少可以拆成多少个完全平方数的和

Input

第一行一个整数TEST，表示数据组数
接下来TEST行，每行一个整数N

Output

TEST行，每行一个正整数，表示N最少可以拆成多少个完全平方数的和。  
0≤N≤10^18 ,TEST<=1000 

Sample Input

5

413357

567226

766291

597007

289611

Sample Output

3

3

3

4

3

拉格朗日四平方和定理：每个正整数均可表示为4个整数的平方和

所以这道题的答案不会超过4

①答案为1：直接对n开更号判断即可

②答案为2：将n进行质因数分解（因为较大所以要用Pollard_rho算法），所有对4取模为3的质因数次数都是偶数

性质：假设n为两个答案为2的数相乘，那么Ans[n]还是等于2

③答案为3：如果不满足其他所有条件，那么答案就为3了hh

④答案为4：n可以表示成4^x(8m+7)(x和m都为整数)的形式

#include<stdio.h>
#include<math.h>
#include<algorithm>
using namespace std;
#define LL long long
int t, cnt, pri[50005], a[50005] = {1,1};
LL fat[101];
LL Multi(LL a, LL b, LL mod)
{
	LL t;
    a %= mod, b %= mod;
    t = a*b-(LL)((double)a/mod*b+0.5)*mod;
	t = (t+mod)%mod;
    return t;
}
/*
LL Multi(LL a, LL b, LL mod)		//和上面一样，但比上面慢很多
{
	LL ans = 0;
	a %= mod;
	while(b)
	{
		if(b%2==1)  ans = (ans+a)%mod, b--;
		else  a = (a+a)%mod, b /= 2;
	}
	return ans;
}
*/
LL Pow(LL a, LL b, LL mod)
{
	LL ans = 1;
	a %= mod;
	while(b)
	{
		if(b&1)  ans = Multi(ans, a, mod), b--;
		else  a = Multi(a, a, mod), b /= 2;
	}
	return ans;
}
LL Gcd(LL a, LL b)
{
	if(b==0)
		return a;
	return Gcd(b, a%b);
}
int Miller_Rabin(LL n)
{
	int i, j, k;
	LL a, x, y, mod;
	if(n==2)  return 1;
	if(n<2 || n%2==0)  return 0;
	k = 0, mod = n-1;
	while(mod%2==0)
	{
		k++;
		mod /= 2;
	}
	for(i=1;i<=10;i++)
	{
		a = rand()%(n-1)+1;
		x = Pow(a, mod, n);
		y = 0;
		for(j=1;j<=k;j++)
		{
			y = Multi(x, x, n);
			if(y==1 && x!=1 && x!=n-1)
				return 0;
			x = y;
		}
		if(y!=1)
			return 0;
	}
	return 1;
}
LL Divi(LL n)
{
	LL i, k, x, y, p, c;
	if(n==1)
		return 1;
	k = 2, p = 1;
	y = x = rand()%n, c = rand()%(n-1)+1;
	for(i=1;p==1;i++)
	{
		x = (Multi(x, x, n)+c)%n;
		p = abs(x-y);
		p = Gcd(n, p);
		if(i==k)
			y = x, k *= 2;
	}
	return p;
}

void Pollard_rho(LL n)
{
	LL p;
	if(n==1)
		return;
	if(Miller_Rabin(n))
		fat[++cnt] = n;
	else
	{
		p = Divi(n);
		Pollard_rho(n/p);
		Pollard_rho(p);
	}
}

int Jud2(LL n)
{
	int i, sum;
	cnt = 0;
	Pollard_rho(n);
	sort(fat+1, fat+1+cnt);
	for(i=1;i<=cnt;i++)
	{
		if(fat[i]%4==3)
		{
			sum = 1;
			while(fat[i]==fat[i+1])
				i++, sum++;
			if(sum%2==1)
				return 0;
		}
	}
	return 1;
}

int main(void)
{
	LL T, n, i, j, t;
	for(i=2;i<50000;i++)
	{
		if(a[i]==0)
			pri[++cnt] = i;
		for(j=i*i;j<=50000;j+=i)
			a[j] = 1;
	}
	scanf("%lld", &T);
	while(T--)
	{
		scanf("%lld", &n);
		t = sqrt(n);
		if(t*t==n)
			printf("1\n");
		else if(Jud2(n))
			printf("2\n");
		else
		{
			while((n&3)==0)
				n /= 4;
			if(n%8==7)
				printf("4\n");
			else
				printf("3\n");
		}
	}
	return 0;
}