poj 3361 Gaussian Prime Factors 高斯素数约数

本文主要是介绍poj 3361 Gaussian Prime Factors 高斯素数约数，希望对大家解决编程问题提供一定的参考价值，需要的开发者们随着小编来一起学习吧！

poj 3361 Gaussian Prime Factors

题意：

在复数 a+bj (a,b为整数)范围内，约数只有 1, -1, a+bj, -(a+bj)的称为高斯素数。求任给正整数N的所有素因子(|b|>a>0)。默认N<2e9

解法：

定理有云：
n≡3(mod 4)者，无因式分解。
n≡1(mod 4)者，可唯一分解为两个共轭高斯整数乘积。

当然做的时候并不知道！
我的思路是，既然可以分解因子，那么可以先分解成素因子，再把素因子分解出来。发现3,7等数是不能分解的，其他的素数(忽略2)都是奇数，要分解成乘积，最后不能有虚部存在，故得分解成共轭高斯整数相乘，即(a+bj)(a-bj)=n，即 a^2+b^2 = n,则a,b一奇一偶。且因为若a+bj不为素数，则可分解为 (c+dj)(e+fj)，那么(a+bj)(a-bj)=(c+dj)(c-dj)(e+fj)(e-fj)=n,此时前两项与后两项都为一个实整数，与n为素数矛盾，故只要分解出来，a+bj与a-bj则一定为n的素因子。
这样就好办了，对所有的N，都先进行素因子分解，然后对每个素因子再分解其高斯素因子。我选择遍历a,b来预处理出所有素数的标准因子a+bj（另一个则为a-bj）。

自己的坑：
低估了45000内素数数量，大约为十分之一个。
使用了set,但是貌似不是这么用的。(我用来验证m-i^2是否为完全平方数，多了一个log(num_pri)的复杂度)。应该是直接检验比较快。

so,学到了STL set:

myset.count(n)  ->得到数量
myset.insert(n) ->插入
myset.find(n)   ->得到迭代器，找不到则是末尾。 .erase(.find(n))

140+行代码，500+ms过。可能是预处理和数据结构比较费时间。看别人的当场暴力分解质因数，暴力分解素数的0ms过。只能说数据太水太水啦！~

我的代码：

#include<iostream>
#include<cstdio>
#include<string>
#include<cstring>
#include<vector>
#include<cmath>
#include<set>
#include<algorithm>
using namespace std;
typedef long long ll;
const int maxn = 45000;
const int maxq = 212;
struct stt{int a,b;
}q[maxn],tmpst;
vector <stt> ans;
set <int> myset;
int prim[5000];
bool isnotpri[maxn];
int prinum = 0;
int prime,sise;
int sq[maxq+1];//squarevoid ini(){memset(prim,0,sizeof prim);memset(isnotpri,false,sizeof isnotpri);memset(q,0,sizeof q);myset.clear();//memset(ans,0,sizeof ans);ans.clear();for(int i = 2; i <= maxq; i++){if(!isnotpri[i]){for(int j = i*i; j < maxn; j += i){isnotpri[j] = 1;}}}for(int i = 2; i <= maxn; i++){if(!isnotpri[i]) prim[prinum++]=i;}q[2].a=1;q[2].b=1;for(int i = 1; i <= maxq; i++){sq[i] = i*i;}int tmp;for(int i = 1; i <= maxq; i++){for(int j = i+1; j <= maxq; j += 2){tmp = sq[i]+sq[j];if(tmp >= maxn) break;if(!isnotpri[tmp]){q[tmp].a = i;q[tmp].b = j;}}}for(int i = 2; i <= maxn; i++){if(q[i].b == 0) q[i].a = i;myset.insert(i*i);}
}
void gener(int m){ans.clear();if(m == 0 || m == 1) return;for(int i = 0; i < prinum; i++){prime = prim[i];if(prime > m) break;if(m % prime == 0){ans.push_back(q[prime]);while(m % prime == 0){m /= prime;}}}if(m != 1){if(m%4==3){tmpst.a = m;tmpst.b = 0;ans.push_back(tmpst);}else{int tmp = (int)sqrt(m+0.1);int tmp2;for(int i = 1; i <= tmp; i++){if(myset.count(m-i*i)!=0){tmpst.a = i;tmpst.b = (int)sqrt(m-i*i+0.1);ans.push_back(tmpst);break;}}}}
}
bool cmp(stt i, stt j){if(i.a == j.a){return i.b < j.b;}return i.a < j.a;
}
void printout(){sise = ans.size();if(sise == 0){printf("\n");return;}sort(ans.begin(),ans.end(),cmp);int aa,bb;for(int i = 0; i < sise; i++){aa = ans[i].a, bb = ans[i].b;if(bb==0){printf(" %d",aa);}else{if(bb==1){printf(" %d+j, %d-j",aa,aa);}else{printf(" %d+%dj, %d-%dj",aa,bb,aa,bb);}}printf("%c",",\n"[i==sise-1]);}
}
int main(){freopen("in.txt","r",stdin);ini();int k = 1;int m;//return 0;while(scanf("%d",&m)!=EOF){printf("Case #%d:",k++);gener(m);printout();}return 0;
}

别人的代码：

from yjCola

#include<cstdio>
#include<cmath>
#include<cstring>
#include<algorithm>
using namespace std;
#define MAXN 1005
struct ele
{int a,b;bool operator < (const ele &c) const{if(a==c.a){if(abs(b)==abs(c.b))return b>c.b;elsereturn abs(b)<abs(c.b);}elsereturn a<c.a;}
}e[MAXN];
int up;
void get(int p)
{if((p-3)%4==0) e[up].a=p,e[up++].b=0;else{for(int i=1;i*i<=p;i++){int j=sqrt(1.0*(p-i*i));if(i*i+j*j==p){e[up].a=i,e[up++].b=j;e[up].a=i,e[up++].b=-j;break;}}}
}
int main()
{int cas=1,n;while(scanf("%d",&n)!=EOF){up=0;for(int i=2;i*i<=n;i++)if(n%i==0){get(i);while(n%i==0) n/=i;}if(n>1) get(n);sort(e,e+up);printf("Case #%d:",cas++);for(int i=0;i<up;i++){printf(" %d",e[i].a);if(e[i].b<0){if(e[i].b==-1) printf("-j");else printf("%dj",e[i].b);}else if(e[i].b>0){if(e[i].b==1) printf("+j");else printf("+%dj",e[i].b);}if(i<up-1) printf(",");}printf("\n");}return 0;
}