什么是最快的分解algorithm？

我已经写了一个试图findAmicable Pairs的程序。 这就要求find合适的数字除数。

这是我目前的sumOfDivisors()方法：

 int sumOfDivisors(int n) { int sum = 1; int bound = (int) sqrt(n); for(int i = 2; i <= 1 + bound; i++) { if (n % i == 0) sum = sum + i + n / i; } return sum; } 

所以我需要做很多因子分解，并且开始成为我应用程序中的真正瓶颈。 我在MAPLE中input了一个很大的数字，并将其快速分解。

什么是更快的分解algorithm之一？

Answer 1

直接从我对这个问题的答案拉。

该方法将工作，但会很慢。 “你的号码有多大？” 决定使用的方法：

• 不到2 ^ 16左右：查询表。
• 不到2 ^ 70左右： Richard Brent修改了Pollard的rhoalgorithm 。
• 小于10 ^ 50： Lenstra椭圆曲线分解
• 小于10 ^ 100： 二次筛
• 10 ^ 100以上： 通用号码筛

Answer 2

标题（和最后一行）中的问题似乎与问题的实际主体毫无关系。 如果你试图find友好的对，或者计算许多数字的除数之和，那么单独分解每个数字（即使采用最快的algorithm）也是一种效率低下的方法。

除数函数 σ(n) = (sum of divisors of n)是一个乘性函数 ：对于相对于素数m和n，我们有σ(mn) = σ(m)σ(n) ，所以

σ（p ₁ ^{k ₁} … p _r ^{k _r} ）= [（p ₁ ^{k ₁ +1} -1）/（p ₁ -1）] … [（p _r ^{k _r +1} -1）/（p _r -1 ）]。

所以你可以使用任何简单的筛子（例如Eratosthenes筛子的扩充版本）来find n数，在这个过程中，所有数字的因式分解都是n。 （例如，当你筛选的时候，存储每个n的最小素因子，然后你可以通过迭代来分解任意数n ）。这将比使用任何单独的分解algorithm几次更快（总体）。

顺便说一句：已经存在几个已知的友好对的列表（例如在这里和MathWorld的链接） – 你也试图扩大logging，或只是为了好玩吗？

Answer 3

Shor的algorithm： http : //en.wikipedia.org/wiki/Shor%27s_algorithm

当然，你需要一台量子计算机：D

Answer 4

我build议从Maple中使用的相同algorithm开始， 二次筛 。

1. select你的奇数n来分解，
2. select一个自然数k ，
3. search所有p <= k ，使得k ^ 2与（n mod p）不一致以获得因子基B = p1，p2，…，pt ，
4. 从r > floor（n）开始search至lesst + 1个值，使得y ^ 2 = r ^ 2 – n在B中都只有因子，
5. 对于刚刚计算的每个y1 ， y2 ，…， y（t + 1） ，生成一个向量v（yi）=（e1，e2，…，et） ，其中ei是通过将模2减去指数pi在yi ，
6. 使用高斯消除来find一些相加的向量给出一个空向量
7. 将x设为前一步find的与yi相关的ri的乘积，并将y设为p1 ^ a * p2 ^ b * p3 ^ c * .. * pt ^ z其中指数是在因子分解中find的指数的一半义
8. 计算d = mcd（xy，n） ，如果1 <d <n，则d是n的非平凡因子，否则从步骤2开始select更大的k。

关于这些algorithm的问题是他们确实暗示了数值计算中的许多理论。

Answer 5

这是Maple中Integer因式分解的一篇文章。

“从一些非常简单的指令开始 – ”在Maple中使整数分解更快“ – 我们在Maple和C的组合中实现了二次筛分解algorithm…”

~pborwein/MITACS/papers/percival.pdf

Answer 6

取决于你的数字有多大。 如果你正在寻找友好的对子，你会做很多因素分解，所以关键可能不是尽快分解，而是要在不同的呼叫之间分享尽可能多的工作。 为了加快审判速度，你可以看一下记忆，和/或预先计算你最关心的人数的平方根。 获得素数因子分解的速度更快，然后计算所有因子的总和，而不是一直循环到每个数字的sqrt（n）。

如果你正在寻找真正的友好对，比如大于2 ^ 64，那么在less数机器上，无论你的因子分解速度有多快，都无法通过分解每一个数字来实现。 你用来find候选人的捷径可能会帮助你分析他们。

Answer 7

1GB内存的更多2015 C ++版本2 ²⁷查询表实现：

 #include <iostream.h> // cerr, cout, and NULL #include <string.h> // memcpy() #define uint unsigned __int32 uint *factors; const uint MAX_F=134217728; // 2^27 void buildFactors(){ factors=new (nothrow) uint [(MAX_F+1)*2]; // 4 * 2 * 2^27 = 2^30 = 1GB if(factors==NULL)return; // not able to allocate enough free memory int i; for(i=0;i<(MAX_F+1)*2;i++)factors[i]=0; //Sieve of Eratosthenese factors[1*2]=1; factors[1*2+1]=1; for(i=2;i*i<=MAX_F;i++){ for(;factors[i*2] && i*i<=MAX_F;i++); factors[i*2]=1; factors[i*2+1]=i; for(int j=2;i*j<=MAX_F;j++){ factors[i*j*2]=i; factors[i*j*2+1]=j; } } for(;i<=MAX_F;i++){ for(;i<=MAX_F && factors[i*2];i++); if(i>MAX_F)return; factors[i*2]=1; factors[i*2+1]=i; } } uint * factor(uint x, int &factorCount){ if(x > MAX_F){factorCount=-1;return NULL;} uint tmp[70], at=x; int i=0; while(factors[at*2]>1){ tmp[i++]=factors[at*2]; cout<<"at:"<<at<<" tmp:"<<tmp[i-1]<<endl; at=factors[at*2+1]; } if(i==0){ cout<<"at:"<<x<<" tmp:1"<<endl; tmp[i++]=1; tmp[i++]=x; }else{ cout<<"at:"<<at<<" tmp:1"<<endl; tmp[i++]=at; } factorCount=i; uint *ret=new (nothrow) uint [factorCount]; if(ret!=NULL) memcpy(ret, tmp, sizeof(uint)*factorCount); return ret; } void main(){ cout<<"Loading factors lookup table"<<endl; buildFactors(); if(factors==NULL){cerr<<"Need 1GB block of free memory"<<endl;return;} int size; uint x=30030; cout<<"\nFactoring: "<<x<<endl; uint *f=factor(x,size); if(size<0){cerr<<x<<" is too big to factor. Choose a number between 1 and "<<MAX_F<<endl;return;} else if(f==NULL){cerr<<"ran out of memory trying to factor "<<x<<endl;return;} cout<<"\nThe factors of: "<<x<<" {"<<f[0]; for(int i=1;i<size;i++) cout<<", "<<f[i]; cout<<"}"<<endl; delete [] f; x=30637; cout<<"\nFactoring: "<<x<<endl; f=factor(x,size); cout<<"\nThe factors of: "<<x<<" {"<<f[0]; for(int i=1;i<size;i++) cout<<", "<<f[i]; cout<<"}"<<endl; delete [] f; delete [] factors; } 

更新：或者在2 ^{月28日之前}牺牲一些简单性

 #include <iostream.h> // cerr, cout, and NULL #include <string.h> // memcpy(), memset() //#define dbg(A) A #ifndef dbg #define dbg(A) #endif #define uint unsigned __int32 #define uint8 unsigned __int8 #define uint16 unsigned __int16 uint * factors; uint8 *factors08; uint16 *factors16; uint *factors32; const uint LIMIT_16 = 514; // First 16-bit factor, 514 = 2*257 const uint LIMIT_32 = 131074;// First 32-bit factor, 131074 = 2*65537 const uint MAX_FACTOR = 268501119; //const uint64 LIMIT_64 = 8,589,934,594; // First 64-bit factor, 2^33+1 const uint TABLE_SIZE = 268435456; // 2^28 => 4 * 2^28 = 2^30 = 1GB 32-bit table const uint o08=1, o16=257 ,o32=65665; //o64=4294934465 // TableSize = 2^37 => 8 * 2^37 = 2^40 1TB 64-bit table // => MaxFactor = 141,733,953,600 /* Layout of factors[] array * Indicies(32-bit) i Value Size AFactorOf(i) * ---------------- ------ ---------- ---------------- * factors[0..128] [1..513] 8-bit factors08[i-o08] * factors[129..65408] [514..131073] 16-bit factors16[i-o16] * factors[65409..268435455] [131074..268501119] 32-bit factors32[i-o32] * * Note: stopping at i*i causes AFactorOf(i) to not always be LargestFactor(i) */ void buildFactors(){ dbg(cout<<"Allocating RAM"<<endl;) factors=new (nothrow) uint [TABLE_SIZE]; // 4 * 2^28 = 2^30 = 1GB if(factors==NULL)return; // not able to allocate enough free memory uint i,j; factors08 = (uint8 *)factors; factors16 = (uint16 *)factors; factors32 = factors; dbg(cout<<"Zeroing RAM"<<endl;) memset(factors,0,sizeof(uint)*TABLE_SIZE); //for(i=0;i<TABLE_SIZE;i++)factors[i]=0; //Sieve of Eratosthenese //8-bit values dbg(cout<<"Setting: 8-Bit Values"<<endl;) factors08[1-o08]=1; for(i=2;i*i<LIMIT_16;i++){ for(;factors08[i-o08] && i*i<LIMIT_16;i++); dbg(cout<<"Filtering: "<<i<<endl;) factors08[i-o08]=1; for(j=2;i*j<LIMIT_16;j++)factors08[i*j-o08]=i; for(;i*j<LIMIT_32;j++)factors16[i*j-o16]=i; for(;i*j<=MAX_FACTOR;j++)factors32[i*j-o32]=i; } for(;i<LIMIT_16;i++){ for(;i<LIMIT_16 && factors08[i-o08];i++); dbg(cout<<"Filtering: "<<i<<endl;) if(i<LIMIT_16){ factors08[i-o08]=1; j=LIMIT_16/i+(LIMIT_16%i>0); for(;i*j<LIMIT_32;j++)factors16[i*j-o16]=i; for(;i*j<=MAX_FACTOR;j++)factors32[i*j-o32]=i; } }i--; dbg(cout<<"Setting: 16-Bit Values"<<endl;) //16-bit values for(;i*i<LIMIT_32;i++){ for(;factors16[i-o16] && i*i<LIMIT_32;i++); factors16[i-o16]=1; for(j=2;i*j<LIMIT_32;j++)factors16[i*j-o16]=i; for(;i*j<=MAX_FACTOR;j++)factors32[i*j-o32]=i; } for(;i<LIMIT_32;i++){ for(;i<LIMIT_32 && factors16[i-o16];i++); if(i<LIMIT_32){ factors16[i-o16]=1; j=LIMIT_32/i+(LIMIT_32%i>0); for(;i*j<=MAX_FACTOR;j++)factors32[i*j-o32]=i; } }i--; dbg(cout<<"Setting: 32-Bit Values"<<endl;) //32-bit values for(;i*i<=MAX_FACTOR;i++){ for(;factors32[i-o32] && i*i<=MAX_FACTOR;i++); factors32[i-o32]=1; for(j=2;i*j<=MAX_FACTOR;j++)factors32[i*j-o32]=i; } for(;i<=MAX_FACTOR;i++){ for(;i<=MAX_FACTOR && factors32[i-o32];i++); if(i>MAX_FACTOR)return; factors32[i-o32]=1; } } uint * factor(uint x, int &factorCount){ if(x > MAX_FACTOR){factorCount=-1;return NULL;} uint tmp[70], at=x; int i=0; while(at>=LIMIT_32 && factors32[at-o32]>1){ tmp[i++]=factors32[at-o32]; dbg(cout<<"at32:"<<at<<" tmp:"<<tmp[i-1]<<endl;) at/=tmp[i-1]; } if(at<LIMIT_32){ while(at>=LIMIT_16 && factors16[at-o16]>1){ tmp[i++]=factors16[at-o16]; dbg(cout<<"at16:"<<at<<" tmp:"<<tmp[i-1]<<endl;) at/=tmp[i-1]; } if(at<LIMIT_16){ while(factors08[at-o08]>1){ tmp[i++]=factors08[at-o08]; dbg(cout<<"at08:"<<at<<" tmp:"<<tmp[i-1]<<endl;) at/=tmp[i-1]; } } } if(i==0){ dbg(cout<<"at:"<<x<<" tmp:1"<<endl;) tmp[i++]=1; tmp[i++]=x; }else{ dbg(cout<<"at:"<<at<<" tmp:1"<<endl;) tmp[i++]=at; } factorCount=i; uint *ret=new (nothrow) uint [factorCount]; if(ret!=NULL) memcpy(ret, tmp, sizeof(uint)*factorCount); return ret; } uint AFactorOf(uint x){ if(x > MAX_FACTOR)return -1; if(x < LIMIT_16) return factors08[x-o08]; if(x < LIMIT_32) return factors16[x-o16]; return factors32[x-o32]; } void main(){ cout<<"Loading factors lookup table"<<endl; buildFactors(); if(factors==NULL){cerr<<"Need 1GB block of free memory"<<endl;return;} int size; uint x=13855127;//25255230;//30030; cout<<"\nFactoring: "<<x<<endl; uint *f=factor(x,size); if(size<0){cerr<<x<<" is too big to factor. Choose a number between 1 and "<<MAX_FACTOR<<endl;return;} else if(f==NULL){cerr<<"ran out of memory trying to factor "<<x<<endl;return;} cout<<"\nThe factors of: "<<x<<" {"<<f[0]; for(int i=1;i<size;i++) cout<<", "<<f[i]; cout<<"}"<<endl; delete [] f; x=30637; cout<<"\nFactoring: "<<x<<endl; f=factor(x,size); cout<<"\nThe factors of: "<<x<<" {"<<f[0]; for(int i=1;i<size;i++) cout<<", "<<f[i]; cout<<"}"<<endl; delete [] f; delete [] factors; }