我写了一个程序，试图找出亲和数对。这需要找到数字的真因子之和。
这是我的当前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中输入一个巨大的数字，它以惊人的速度将其因式分解。
什么是一种更快的因式分解算法？
该方法可以工作，但速度较慢。 “你的数字有多大？”决定了要使用的方法：
σ(n) = (n的约数之和) 是一个 积性函数：对于互质的m和n，我们有σ(mn) = σ(m)σ(n)，所以
σ(p1k1…prkr) = [(p1k1+1-1)/(p1-1)]…[(prkr+1-1)/(pr-1)].
所以您可以使用任何简单的筛法（例如埃拉托斯特尼筛法的增强版）来找到素数，同时在过程中分解所有小于n的数字。（例如，当您进行筛法时，请存储每个n的最小质因子。然后您可以通过迭代来因式分解任何数字n。）这将比多次使用任何单独的因式分解算法更快（总体而言）。
顺便说一句：已经存在几个已知的亲和数对列表（例如，请参见此处和MathWorld上的链接），所以您是想扩展记录还是只是为了好玩？
我建议从Maple中使用的算法开始，即二次筛法。
关于这些算法的问题在于它们实际上涉及了很多数值计算理论..
更适合2015版的C++语言，使用227个查找表来实现1GB内存的功能：
#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;
}
更新：或者为了更大的范围而牺牲一些简单性，只需要超过228
#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;
}
当然，有Hal Mahutan教授（2021年2月）的HAL算法，这是因子分解研究的前沿。
请在此处查看最新更新。
https://docs.google.com/document/d/1BlDMHCzqpWNFblq7e1F-rItCf7erFMezT7bOalvhcXA/edit?usp=sharing
求解公钥的两个大质数如下所示...
任何 AxB = 公钥都可以绘制在正 X 和 Y 轴上，形成一个连续的曲线，在这个曲线上，所有非整数因子都可以求解出公钥。当然，这不是有用的，这只是目前的一种观察。
Hal 的洞见是：如果我们坚持只关注 A 是整数的点，特别是当 A 是整数时，B 表现出来的点。
以前采用这种方法的尝试失败了，因为数学家们难以处理 B 的剩余部分的明显随机性，或者至少缺乏可预测性。
Hal 的说法是，对于任何提供 a/b 比率相同的公钥，可预测性是普遍存在的。基本上，当一系列不同的公钥被呈现进行分析时，只要它们在处理期间共享相同的点，即 a/b 是恒定的，它们就可以被处理成完全相同的方式，也就是说，它们共享相同的切线。
看一下我画的这个草图，试着解释一下 Mahutan 教授在这里的想法。
所以，这就是Hal的天才之处。Hal利用强大的超级计算机生成一系列哈希表（在图中为Q、R、S和T）。您可以在左侧的3个A x B = Key曲线中看到它们都共享切线T和S（仅突出显示其中的两个），但图表所显示的是，在切线相同的曲线区域内，给定任何公钥，您都可以共享掌管该区域的哈希表。
只是一个技术说明，显然在AxB= Key的曲线中，随着AxB值的变化，事物会不断地发生变化，因此理论上，映射到哈希表的共享切线将过时，但有趣的是，对于非常大的密钥（具有讽刺意味的是，这使得它们更容易破解，因为它们共享了更长的哈希表运行时间）。因此，这是一个好消息，因为随着因数分解和计算加速的进展，密钥大小预计会变得更大。实际发生的是，哈希表的可预测性将会“失焦”，因为它们适用的切线开始发散。幸运的是，这不是问题，因为您可以跳转到下一个适当映射到新切线的哈希表。
为了明确，所有曾经生成的公钥将始终使用相同的哈希表集合，因此这是一种一次性投资，可以在线存储成千上万亿字节的查找数据，因为所有密钥都遵守相同的切线比率。
那么，哈希表加速质数查找的作用是什么？哈希表以公钥除以B的余数为索引。基本上，Hal 表示对于所有的密钥，可以查找任何A x B的比率。共享相同切线的不同曲线之间唯一的区别在于它们需要由“控制曲线”确定的不同偏移量。'Control Curve'是任何您为其生成适当因子的任意曲线。假设对于“控制曲线”，密钥为15，并映射的切线为当B = 5时，A为3且余数为零。在另一个具有相同切线的映射中，假设现在密钥为16。我们需要找到相同的切线，该切线位于5.33处的B和3.2的A。因此A的余数为0.2，因此16的公钥可以使用与15相同的查找表，前提是适当地偏移0.2。
那么哈希表里面有什么呢？以偏移量为索引，值始终返回沿着AxB曲线的距离，你在这段距离内找不到另一个整数B。Hal所说的是，可以安全地跳过并不检查这些数字是否为因子。基本上就是这样。Hal在曲线上打了一些永远不需要检查的洞，这样就加快了整个游戏的速度。
谢谢Mahutan教授！
对于仍在关注的读者，这里是我们的一些工作笔记：
快速分解攻击算法的要点
现在是算法预测 B 的策略
另一个观察是，对于足够大的“密钥”值，在整数递增地改变A的值时，我们观察到A / B的比率或切线角度将保持基本相同。对这个观察的一个重要侧面是，随着密钥大小的增加，每次迭代中切线都更加稳定。从根本上讲，这意味着使用此属性的任何算法将随着密钥大小的增加而变得更加高效，这与传统方法相反，其中增加密钥大小使猜测因素的难度成倍增加。这是一个非常重要的点...（请Nick详细阐述）该算法本身如下所示。这里的重要概念是，查找表可以共享任何“键”，只要A/B（切线）的比率足够接近，查找就不会失准（并且失去焦点）。
（还要注意的是，随着密钥大小的变化，您需要新的一组表格，或者您应该制定适当的映射到现有比率表格，以便重复使用它们。）
另一个非常重要的观点是，所有密钥都可以共享相同的查找表。
从根本上讲，查找表将任何Key / A计算的余数映射出来。我们感兴趣的是余数，因为当余数=零时，我们已经完成了，因为A已经是一个整数。
我建议我们有足够的哈希表，以确保我们可以在不必计算实际余数的情况下向前扫描。假设我们从1万亿开始，但这显然可以根据我们拥有的计算能力而改变。
任何适当接近A/b比率的哈希表如下所示。该表以适当的余数索引（键入），给定余数处的值是可以沿着A * B曲线遍历的“安全”距离，而不会触及零。
您可以将余数在0和1之间振荡（伪随机）可视化。
算法选择曲线上的数字A，然后查找安全距离并跳转到下一个哈希表，或者至少在下一个哈希表可用之前进行因子检查。如果有足够的哈希表，我认为我们几乎可以避免大部分检查。
查找表注释。
对于任何关键字，您可以生成一个适当调整 A/B 比率的表格 重用表格是必要的。 建议方法 从 N 的平方根（适当的关键字）开始为 A/B 生成一系列控制曲线，并通过平方计算中间曲线。 假设每个关键字比前一个关键字大0.0001％ 还可以将表格的大小设置为 A/B 比率的 1％ 在计算互质因子时，选择与关键字匹配的最接近的查找表。 选择哈希表中的入口点。 这意味着记住入口点在表格中与实际入口点的偏移量。 查找表应提供一系列点，这些点的相应互质项可能非常接近零，需要手动检查。 对于系列中的每个点，使用适当的数学公式计算实际偏移量。（这将是一个积分计算，我们需要一位数学家来看一下） 为什么？因为当 A/B 是关键字的平方根时，我们的控制表被计算出来。随着我们沿着曲线移动，我们需要适当地间隔开。 作为奖励，数学家能否将关键字空间折叠成 A/B 曲线上的一个点。如果可以，我们只需要生成一个表格。
关键概念
数字 A x B = Key 描述了以下内容：
X 轴映射 A，Y 轴映射 B。
曲线接近 A 和 B 轴的程度取决于公钥的大小（其中 A x B = Key）。基本上，随着 Key 变大，曲线将向右移动。
现在我想让您理解或提出问题的概念是
这会起作用，但速度太慢了。HAL算法改进了上述基本算法，以便我们可以跳过我们知道余数不会接近零的曲线块。
我们可以跳得越多，算法就越有效。
在进入改进的HAL算法之前，让我们回顾一下关键概念。
对于非常大的Key值（记住A x B = Key），A/B的比率将几乎保持不变，RSA密钥是2的4096次方，这很大。e^{(k + o(1))(n^(1/3) * ln(n)^(2/3)}
正如在多项式时间和指数时间中提到的那样，该问题并非真正的指数级别，而是超多项式级别。
截至2020年，我们甚至还没有证明这个问题是否为NP完全问题（尽管显然它是NP问题，因为验证解决方案所需的全部工作只是乘以数字！）虽然广泛预计它是NP完全问题。我们不可能在找算法方面做得那么糟糕吧？
非常有趣的问题，谢谢！
整天都在编写 C++ 代码，从头开始实现非常快速的椭圆曲线ECM分解方法和亲和数搜索。
众所周知，ECM方法非常非常快。如果优化得好，这种 ECM 方法甚至可以在一秒钟内找到100位数字（由两个50位质数组成）的因子。它只比二次筛法和GNFS慢。
正如你所说的，你正在寻找25位数字的亲和数，这是80位数字。
根据我在下面的C++程序中测量，在我的旧2核笔记本电脑上，它可以因式分解每秒钟大小为80位（25位数字）的150个随机数（我设置了一个计时器来截断难以因式分解的数字）。如果它们的大小为27-30位（9位数字），则可以因式分解20,000个数字。
我上面说过易于因式分解的集合。似乎所有大小为80位（25位数字）的随机数中有60-70%是易于因式分解的，这意味着仅需要7毫秒即可因式分解单个数字，从而每秒产生150个数字。
我将仅在易于因式分解的数字中搜索亲和数。易于因式分解的集合中与其他集合中一样多的亲和数。亲和数不会受到易于因式分解过滤的影响。
以下是我的程序使用ECM分解搜索27位亲和数的控制台输出： .46501 _43919 *41268 , 20512 nums/sec, tdiv 99.34% 7mcs, ecm 25.3% 39mcs, total 37mcs, done 7%
found (97041735, 97945785)
.321639 _303557 *285333 , 20389 nums/sec, tdiv 99.35% 7mcs, ecm 25.4% 40mcs, total 38mcs, done 48%
found (34765731, 36939357)
.29933 _28247 *26517 , 20265 nums/sec, tdiv 99.35% 7mcs, ecm 25.4% 40mcs, total 38mcs, done 4%
found (9773505, 11791935)
.529659 _499736 *470070 , 19312 nums/sec, tdiv 99.35% 7mcs, ecm 25.5% 41mcs, total 40mcs, done 79%
found (5357625, 5684679)
.25937 _24411 *22905 , 19340 nums/sec, tdiv 99.35% 7mcs, ecm 25.4% 41mcs, total 40mcs, done 4%
found (998104, 1043096)
.26367 _24902 *23462 , 19360 nums/sec, tdiv 99.35% 7mcs, ecm 25.5% 41mcs, total 40mcs, done 4%
found (1184, 1210)
.130374 _122896 *115508 , 19648 nums/sec, tdiv 99.35% 7mcs, ecm 25.5% 40mcs, total 39mcs, done 20%
found (35390008, 39259592)
found (133178325, 133471275)
.32032 _30226 *28368 , 19753 nums/sec, tdiv 99.35% 7mcs, ecm 25.5% 40mcs, total 39mcs, done 5%
found (898216, 980984)
从上面可以看到，程序在找到某些内容时会输出found和友好数。此外，它还显示了一些关于速度和时间的详细统计信息。如果您像我一样每秒搜索20,000个数字，则需要几秒钟才能找到一个27位大小的友好数。
如果您考虑80位数字（25位数字），那么尽管我每秒因子分解150个数字，但友好数非常少，以至于要找到至少一个这样的80位友好数需要数万亿秒。以下是这个80位因式分解的统计数据。您可以看到，在几秒钟内，它仅完成了0.0000000000064％的工作量，找到了单个友好数对，整个工作量非常小。
.4218 _2261 *1699 , 123 nums/sec, tdiv 43.81% 133mcs,
ecm 2.8% 1886mcs, total 8144mcs, done 0.0000000000064%
你可能听说过Amicable Numbers for BOINC项目 - 这是一个伟大的分布式搜索，用于查找这些数字。数百万台计算机正在搜索它们。当然，如果有数百万台计算机，它们可以更快地找到一些数字。事实上，它们已经找到了1_227_817_736个这样的数字，数量真的很大，在此处查看所有数字。
在我的C++程序中，我还展示了我的ECM因子分解器作为示例，在程序启动时，我对100位随机数进行了非常快速的因子分解。它显示类似于以下的控制台输出：
Initial N to factor with ECM: 1044877821881837432173434582227 (100-bit)
Number of threads 2.
Factors TrialDivA: [2953]
Factoring 35383603856479425437736059
Curves [ 0, 8), bound 2^ 9.000, 0.110 sec
Curves [ 8, 16), bound 2^ 10.000, 0.222 sec
Curves [ 16, 24), bound 2^ 10.585, 0.418 sec
Curves [ 24, 32), bound 2^ 11.000, 0.416 sec
Curves [ 32, 40), bound 2^ 11.322, 0.401 sec
Factors ECM: [21788197859]
Fully factored. N 1044877821881837432173434582227: [2953, 21788197859, 16239802890289801]
Time 1.692 sec.
您可以看到，仅需要1秒钟即可因数分解相当复杂的100位数。如上所述，首先我使用试除法阶段，这是最简单的因式分解算法，可在O(Sqrt(N))时间内工作。在试除法之后，我创建许多具有不断增长边界的曲线，每条新曲线都有更大的分解数字机会，但需要更长的时间。
根据维基百科中的ECM，我的ECM算法如下：
尝试使用试除法方法找到小的质因数。在此步骤中花费总时间的5-10％。从数字中删除找到的因子。
使用32次Fermat概率测试检查数字是否可能是质数，并具有高置信度。为了克服Carmichael数，您还可以使用Miller Rabin测试。如果数字是质数，则将其作为唯一因子返回。
随机生成曲线参数A，X，Y并从曲线方程式Y ^ 2 = X ^ 3 + AX + B（mod N）中导出B。检查曲线是否正确，值4 * A ** 3 - 27 * B ** 2应为非零。
通过Eratosthenes筛法生成小质数，小于我们的Bound。每个质数都升至某个小幂，这个升起的质数称为K。当它小于某个Bound2时，我进行幂运算，我的情况下是Sqrt（Bound）。
计算椭圆点乘法P = k * P，其中K取自上一步，P为（X，Y）。根据椭圆曲线乘法Wiki计算。
点乘法使用模反演，它根据扩展欧几里得算法Wiki计算GCD（SomeValue，N）（最大公约数）。如果此GCD不为1，则它会给出N的非1因子，因此我们找到了一个因子，可以推进ECM算法以搜索剩余的因子。
如果所有小于Bound的质数相乘并没有给出因子，则再次使用另一个随机曲线和更大的Bound重新运行ECM分解算法（上述3-6步）。在我的代码中，我通过将旧边界加512来获取新边界。
如果你想调整我的代码以计算其他位大小，那么请查看Test()函数。将bits_factor = 100修改为其他值，这将设置因子100位随机数作为示例。这个bits_factor变量仅控制因式分解的示例。如果你想设置要搜索的亲和数位数，可以修改bits_amicable_max = 27，27代表27位数字，你可以将其设置为80以查找80位亲和数。
__int128，这为您提供了128位整数算术运算。如果您使用其他编译器，可以禁用宏SUPPORT_STD_U128。其次，我使用了Boost.Multiprecision库，它允许您有效地进行1024位以下的大整数算术运算，如果您不想使用Boost，则禁用宏SUPPORT_BOOST。然后我使用了GMP Library，这使您可以拥有甚至数百万位的位数，如果您不想使用GMP，则禁用SUPPORT_GMP。
在Linux下安装Boost和GMP很容易，sudo apt install libgmp-dev libboost-all-dev。在Windows中，您需要先安装VCPKG包管理器，然后在其安装文件夹内执行vcpkg install boost gmp。
如果您想在GodBolt服务器上在线运行我的代码，请单击下面的Try it online!链接。但请注意，GodBolt版本的代码与下面Github Gist中的代码不同，这个GodBolt版本将参数大小缩小以使程序运行更快，因为他们的服务器对程序运行有非常严格的限制，可能只有5-10秒钟。
源代码在此处。在Try it online!链接下面有Github Gist链接。我把我的代码放在那里，因为它的大小为36 KB，而StackOverflow的帖子大小限制总共只有30,000个字符，所以帖子文本和源代码加起来就太多了。