将数字分解为 2 个质数余因子
Decompose a number into 2 prime co-factors
Telegram Authentication 的要求之一是将给定数字分解为 2 个素数余因子。特别是 P*Q = N, where N < 2^63
我们怎样才能找到更小的素余因子，使得P < square_root(N)
我的建议：
1) 预先计算从 3 到 2^31.5 的素数，然后测试是否 N mod P = 0
2) 找到一个算法来测试素数（但是我们还是要测试N mod P =0）
是否有适合这种情况的素数算法？
呃！我刚把这个程序放进去，然后意识到你已经用 C# 标记了你的问题。这是 C++，我几年前写的 Pollard Rho 的一个版本，并在 SO 上发布以帮助其他人理解它。分解半素数比试除法快很多倍。正如我所说，我很遗憾它是 C++ 而不是 C#，但是您应该能够理解这个概念，甚至可以很容易地移植它。作为奖励，.NET 库有一个用于处理任意大整数的命名空间，我的 C++ 实现要求我为它们寻找第三方库。无论如何，即使在 C# 中，下面的程序也会在不到 1 秒的时间内将 2^63 阶半素数分解为 2 个素数。甚至还有比这更快的算法，但它们要复杂得多。
#include <string>
#include <stdio.h>
#include <iostream>
#include "BigIntegerLibrary.hh"
typedef BigInteger BI;
typedef BigUnsigned BU;
using std::string;
using std::cin;
using std::cout;
BU pollard(BU &numberToFactor);
BU gcda(BU differenceBetweenCongruentFunctions, BU numberToFactor);
BU f(BU &x, BU &numberToFactor, int &increment);
void initializeArrays();
BU getNumberToFactor ();
void factorComposites();
bool testForComposite (BU &num);
BU primeFactors[1000];
BU compositeFactors[1000];
BU tempFactors [1000];
int primeIndex;
int compositeIndex;
int tempIndex;
int numberOfCompositeFactors;
bool allJTestsShowComposite;
int main ()
{
while(1)
{
primeIndex=0;
compositeIndex=0;
tempIndex=0;
initializeArrays();
compositeFactors[0] = getNumberToFactor();
cout<<"\n\n";
if (compositeFactors[0] == 0) return 0;
numberOfCompositeFactors = 1;
factorComposites();
}
}
void initializeArrays()
{
for (int i = 0; i<1000;i++)
{
primeFactors[i] = 0;
compositeFactors[i]=0;
tempFactors[i]=0;
}
}
BU getNumberToFactor ()
{
std::string s;
std::cout<<"Enter the number for which you want a prime factor, or 0 to quit: ";
std::cin>>s;
return stringToBigUnsigned(s);
}
void factorComposites()
{
while (numberOfCompositeFactors!=0)
{
compositeIndex = 0;
tempIndex = 0;
// This while loop finds non-zero values in compositeFactors.
// If they are composite, it factors them and puts one factor in tempFactors,
// then divides the element in compositeFactors by the same amount.
// If the element is prime, it moves it into tempFactors (zeros the element in compositeFactors)
while (compositeIndex < 1000)
{
if(compositeFactors[compositeIndex] == 0)
{
compositeIndex++;
continue;
}
if(testForComposite(compositeFactors[compositeIndex]) == false)
{
tempFactors[tempIndex] = compositeFactors[compositeIndex];
compositeFactors[compositeIndex] = 0;
tempIndex++;
compositeIndex++;
}
else
{
tempFactors[tempIndex] = pollard (compositeFactors[compositeIndex]);
compositeFactors[compositeIndex] /= tempFactors[tempIndex];
tempIndex++;
compositeIndex++;
}
}
compositeIndex = 0;
// This while loop moves all remaining non-zero values from compositeFactors into tempFactors
// When it is done, compositeFactors should be all 0 value elements
while (compositeIndex < 1000)
{
if (compositeFactors[compositeIndex] != 0)
{
tempFactors[tempIndex] = compositeFactors[compositeIndex];
compositeFactors[compositeIndex] = 0;
tempIndex++;
compositeIndex++;
}
else compositeIndex++;
}
compositeIndex = 0;
tempIndex = 0;
// This while loop checks all non-zero elements in tempIndex.
// Those that are prime are shown on screen and moved to primeFactors
// Those that are composite are moved to compositeFactors
// When this is done, all elements in tempFactors should be 0
while (tempIndex<1000)
{
if(tempFactors[tempIndex] == 0)
{
tempIndex++;
continue;
}
if(testForComposite(tempFactors[tempIndex]) == false)
{
primeFactors[primeIndex] = tempFactors[tempIndex];
cout<<primeFactors[primeIndex]<<"\n";
tempFactors[tempIndex]=0;
primeIndex++;
tempIndex++;
}
else
{
compositeFactors[compositeIndex] = tempFactors[tempIndex];
tempFactors[tempIndex]=0;
compositeIndex++;
tempIndex++;
}
}
compositeIndex=0;
numberOfCompositeFactors=0;
// This while loop just checks to be sure there are still one or more composite factors.
// As long as there are, the outer while loop will repeat
while(compositeIndex<1000)
{
if(compositeFactors[compositeIndex]!=0) numberOfCompositeFactors++;
compositeIndex ++;
}
}
return;
}
// The following method uses the Miller-Rabin primality test to prove with 100% confidence a given number is composite,
// or to establish with a high level of confidence -- but not 100% -- that it is prime
bool testForComposite (BU &num)
{
BU confidenceFactor = 101;
if (confidenceFactor >= num) confidenceFactor = num-1;
BU a,d,s, nMinusOne;
nMinusOne=num-1;
d=nMinusOne;
s=0;
while(modexp(d,1,2)==0)
{
d /= 2;
s++;
}
allJTestsShowComposite = true; // assume composite here until we can prove otherwise
for (BI i = 2 ; i<=confidenceFactor;i++)
{
if (modexp(i,d,num) == 1)
continue; // if this modulus is 1, then we cannot prove that num is composite with this value of i, so continue
if (modexp(i,d,num) == nMinusOne)
{
allJTestsShowComposite = false;
continue;
}
BU exponent(1);
for (BU j(0); j.toInt()<=s.toInt()-1;j++)
{
exponent *= 2;
if (modexp(i,exponent*d,num) == nMinusOne)
{
// if the modulus is not right for even a single j, then break and increment i.
allJTestsShowComposite = false;
continue;
}
}
if (allJTestsShowComposite == true) return true; // proven composite with 100% certainty, no need to continue testing
}
return false;
/* not proven composite in any test, so assume prime with a possibility of error =
(1/4)^(number of different values of i tested). This will be equal to the value of the
confidenceFactor variable, and the "witnesses" to the primality of the number being tested will be all integers from
2 through the value of confidenceFactor.
Note that this makes this primality test cryptographically less secure than it could be. It is theoretically possible,
if difficult, for a malicious party to pass a known composite number for which all of the lowest n integers fail to
detect that it is composite. A safer way is to generate random integers in the outer "for" loop and use those in place of
the variable i. Better still if those random numbers are checked to ensure no duplicates are generated.
*/
}
BU pollard(BU &n)
{
if (n == 4) return 2;
BU x = 2;
BU y = 2;
BU d = 1;
int increment = 1;
while(d==1||d==n||d==0)
{
x = f(x,n, increment);
y = f(y,n, increment);
y = f(y,n, increment);
if (y>x)
{
d = gcda(y-x, n);
}
else
{
d = gcda(x-y, n);
}
if (d==0)
{
x = 2;
y = 2;
d = 1;
increment++; // This changes the pseudorandom function we use to increment x and y
}
}
return d;
}
BU gcda(BU a, BU b)
{
if (a==b||a==0)
return 0; // If x==y or if the absolute value of (x-y) == the number to be factored, then we have failed to find
// a factor. I think this is not proof of primality, so the process could be repeated with a new function.
// For example, by replacing x*x+1 with x*x+2, and so on. If many such functions fail, primality is likely.
BU currentGCD = 1;
while (currentGCD!=0) // This while loop is based on Euclid's algorithm
{
currentGCD = b % a;
b=a;
a=currentGCD;
}
return b;
}
BU f(BU &x, BU &n, int &increment)
{
return (x * x + increment) % n;
}
Pollard 的 Rho 算法 [VB.Net]
查找 P 非常快，其中 P*Q = N，对于 N < 2^63
Dim rnd As New System.Random
Function PollardRho(n As BigInteger) As BigInteger
If n Mod 2 = 0 Then Return 2
Dim x As BigInteger = rnd.Next(1, 1000)
Dim c As BigInteger = rnd.Next(1, 1000)
Dim g As BigInteger = 1
Dim y = x
While g = 1
x = ((x * x) Mod n + c) Mod n
y = ((y * y) Mod n + c) Mod n
y = ((y * y) Mod n + c) Mod n
g = gcd(BigInteger.Abs(x - y), n)
End While
Return g
End Function
Function gcd(a As BigInteger, b As BigInteger) As BigInteger
Dim r As BigInteger
While b <> 0
r = a Mod b
a = b
b = r
End While
Return a
End Function
Richard Brent 的算法 [VB.Net] 这甚至更快。
Function Brent(n As BigInteger) As BigInteger
If n Mod 2 = 0 Then Return 2
Dim y As BigInteger = rnd.Next(1, 1000)
Dim c As BigInteger = rnd.Next(1, 1000)
Dim m As BigInteger = rnd.Next(1, 1000)
Dim g As BigInteger = 1
Dim r As BigInteger = 1
Dim q As BigInteger = 1
Dim x As BigInteger = 0
Dim ys As BigInteger = 0
While g = 1
x = y
For i = 1 To r
y = ((y * y) Mod n + c) Mod n
Next
Dim k = New BigInteger(0)
While (k < r And g = 1)
ys = y
For i = 1 To BigInteger.Min(m, r - k)
y = ((y * y) Mod n + c) Mod n
q = q * (BigInteger.Abs(x - y)) Mod n
Next
g = gcd(q, n)
k = k + m
End While
r = r * 2
End While
If g = n Then
While True
ys = ((ys * ys) Mod n + c) Mod n
g = gcd(BigInteger.Abs(x - ys), n)
If g > 1 Then
Exit While
End If
End While
End If
Return g
End Function
Telegram Authentication 的要求之一是将给定数字分解为 2 个素数余因子。特别是 P*Q = N, where N < 2^63
我们怎样才能找到更小的素余因子，使得P < square_root(N)
我的建议：
1) 预先计算从 3 到 2^31.5 的素数，然后测试是否 N mod P = 0
2) 找到一个算法来测试素数（但是我们还是要测试N mod P =0）
是否有适合这种情况的素数算法？
呃！我刚把这个程序放进去，然后意识到你已经用 C# 标记了你的问题。这是 C++，我几年前写的 Pollard Rho 的一个版本，并在 SO 上发布以帮助其他人理解它。分解半素数比试除法快很多倍。正如我所说，我很遗憾它是 C++ 而不是 C#，但是您应该能够理解这个概念，甚至可以很容易地移植它。作为奖励，.NET 库有一个用于处理任意大整数的命名空间，我的 C++ 实现要求我为它们寻找第三方库。无论如何，即使在 C# 中，下面的程序也会在不到 1 秒的时间内将 2^63 阶半素数分解为 2 个素数。甚至还有比这更快的算法，但它们要复杂得多。
#include <string>
#include <stdio.h>
#include <iostream>
#include "BigIntegerLibrary.hh"
typedef BigInteger BI;
typedef BigUnsigned BU;
using std::string;
using std::cin;
using std::cout;
BU pollard(BU &numberToFactor);
BU gcda(BU differenceBetweenCongruentFunctions, BU numberToFactor);
BU f(BU &x, BU &numberToFactor, int &increment);
void initializeArrays();
BU getNumberToFactor ();
void factorComposites();
bool testForComposite (BU &num);
BU primeFactors[1000];
BU compositeFactors[1000];
BU tempFactors [1000];
int primeIndex;
int compositeIndex;
int tempIndex;
int numberOfCompositeFactors;
bool allJTestsShowComposite;
int main ()
{
while(1)
{
primeIndex=0;
compositeIndex=0;
tempIndex=0;
initializeArrays();
compositeFactors[0] = getNumberToFactor();
cout<<"\n\n";
if (compositeFactors[0] == 0) return 0;
numberOfCompositeFactors = 1;
factorComposites();
}
}
void initializeArrays()
{
for (int i = 0; i<1000;i++)
{
primeFactors[i] = 0;
compositeFactors[i]=0;
tempFactors[i]=0;
}
}
BU getNumberToFactor ()
{
std::string s;
std::cout<<"Enter the number for which you want a prime factor, or 0 to quit: ";
std::cin>>s;
return stringToBigUnsigned(s);
}
void factorComposites()
{
while (numberOfCompositeFactors!=0)
{
compositeIndex = 0;
tempIndex = 0;
// This while loop finds non-zero values in compositeFactors.
// If they are composite, it factors them and puts one factor in tempFactors,
// then divides the element in compositeFactors by the same amount.
// If the element is prime, it moves it into tempFactors (zeros the element in compositeFactors)
while (compositeIndex < 1000)
{
if(compositeFactors[compositeIndex] == 0)
{
compositeIndex++;
continue;
}
if(testForComposite(compositeFactors[compositeIndex]) == false)
{
tempFactors[tempIndex] = compositeFactors[compositeIndex];
compositeFactors[compositeIndex] = 0;
tempIndex++;
compositeIndex++;
}
else
{
tempFactors[tempIndex] = pollard (compositeFactors[compositeIndex]);
compositeFactors[compositeIndex] /= tempFactors[tempIndex];
tempIndex++;
compositeIndex++;
}
}
compositeIndex = 0;
// This while loop moves all remaining non-zero values from compositeFactors into tempFactors
// When it is done, compositeFactors should be all 0 value elements
while (compositeIndex < 1000)
{
if (compositeFactors[compositeIndex] != 0)
{
tempFactors[tempIndex] = compositeFactors[compositeIndex];
compositeFactors[compositeIndex] = 0;
tempIndex++;
compositeIndex++;
}
else compositeIndex++;
}
compositeIndex = 0;
tempIndex = 0;
// This while loop checks all non-zero elements in tempIndex.
// Those that are prime are shown on screen and moved to primeFactors
// Those that are composite are moved to compositeFactors
// When this is done, all elements in tempFactors should be 0
while (tempIndex<1000)
{
if(tempFactors[tempIndex] == 0)
{
tempIndex++;
continue;
}
if(testForComposite(tempFactors[tempIndex]) == false)
{
primeFactors[primeIndex] = tempFactors[tempIndex];
cout<<primeFactors[primeIndex]<<"\n";
tempFactors[tempIndex]=0;
primeIndex++;
tempIndex++;
}
else
{
compositeFactors[compositeIndex] = tempFactors[tempIndex];
tempFactors[tempIndex]=0;
compositeIndex++;
tempIndex++;
}
}
compositeIndex=0;
numberOfCompositeFactors=0;
// This while loop just checks to be sure there are still one or more composite factors.
// As long as there are, the outer while loop will repeat
while(compositeIndex<1000)
{
if(compositeFactors[compositeIndex]!=0) numberOfCompositeFactors++;
compositeIndex ++;
}
}
return;
}
// The following method uses the Miller-Rabin primality test to prove with 100% confidence a given number is composite,
// or to establish with a high level of confidence -- but not 100% -- that it is prime
bool testForComposite (BU &num)
{
BU confidenceFactor = 101;
if (confidenceFactor >= num) confidenceFactor = num-1;
BU a,d,s, nMinusOne;
nMinusOne=num-1;
d=nMinusOne;
s=0;
while(modexp(d,1,2)==0)
{
d /= 2;
s++;
}
allJTestsShowComposite = true; // assume composite here until we can prove otherwise
for (BI i = 2 ; i<=confidenceFactor;i++)
{
if (modexp(i,d,num) == 1)
continue; // if this modulus is 1, then we cannot prove that num is composite with this value of i, so continue
if (modexp(i,d,num) == nMinusOne)
{
allJTestsShowComposite = false;
continue;
}
BU exponent(1);
for (BU j(0); j.toInt()<=s.toInt()-1;j++)
{
exponent *= 2;
if (modexp(i,exponent*d,num) == nMinusOne)
{
// if the modulus is not right for even a single j, then break and increment i.
allJTestsShowComposite = false;
continue;
}
}
if (allJTestsShowComposite == true) return true; // proven composite with 100% certainty, no need to continue testing
}
return false;
/* not proven composite in any test, so assume prime with a possibility of error =
(1/4)^(number of different values of i tested). This will be equal to the value of the
confidenceFactor variable, and the "witnesses" to the primality of the number being tested will be all integers from
2 through the value of confidenceFactor.
Note that this makes this primality test cryptographically less secure than it could be. It is theoretically possible,
if difficult, for a malicious party to pass a known composite number for which all of the lowest n integers fail to
detect that it is composite. A safer way is to generate random integers in the outer "for" loop and use those in place of
the variable i. Better still if those random numbers are checked to ensure no duplicates are generated.
*/
}
BU pollard(BU &n)
{
if (n == 4) return 2;
BU x = 2;
BU y = 2;
BU d = 1;
int increment = 1;
while(d==1||d==n||d==0)
{
x = f(x,n, increment);
y = f(y,n, increment);
y = f(y,n, increment);
if (y>x)
{
d = gcda(y-x, n);
}
else
{
d = gcda(x-y, n);
}
if (d==0)
{
x = 2;
y = 2;
d = 1;
increment++; // This changes the pseudorandom function we use to increment x and y
}
}
return d;
}
BU gcda(BU a, BU b)
{
if (a==b||a==0)
return 0; // If x==y or if the absolute value of (x-y) == the number to be factored, then we have failed to find
// a factor. I think this is not proof of primality, so the process could be repeated with a new function.
// For example, by replacing x*x+1 with x*x+2, and so on. If many such functions fail, primality is likely.
BU currentGCD = 1;
while (currentGCD!=0) // This while loop is based on Euclid's algorithm
{
currentGCD = b % a;
b=a;
a=currentGCD;
}
return b;
}
BU f(BU &x, BU &n, int &increment)
{
return (x * x + increment) % n;
}
Pollard 的 Rho 算法 [VB.Net]
查找 P 非常快，其中 P*Q = N，对于 N < 2^63
Dim rnd As New System.Random
Function PollardRho(n As BigInteger) As BigInteger
If n Mod 2 = 0 Then Return 2
Dim x As BigInteger = rnd.Next(1, 1000)
Dim c As BigInteger = rnd.Next(1, 1000)
Dim g As BigInteger = 1
Dim y = x
While g = 1
x = ((x * x) Mod n + c) Mod n
y = ((y * y) Mod n + c) Mod n
y = ((y * y) Mod n + c) Mod n
g = gcd(BigInteger.Abs(x - y), n)
End While
Return g
End Function
Function gcd(a As BigInteger, b As BigInteger) As BigInteger
Dim r As BigInteger
While b <> 0
r = a Mod b
a = b
b = r
End While
Return a
End Function
Richard Brent 的算法 [VB.Net] 这甚至更快。
Function Brent(n As BigInteger) As BigInteger
If n Mod 2 = 0 Then Return 2
Dim y As BigInteger = rnd.Next(1, 1000)
Dim c As BigInteger = rnd.Next(1, 1000)
Dim m As BigInteger = rnd.Next(1, 1000)
Dim g As BigInteger = 1
Dim r As BigInteger = 1
Dim q As BigInteger = 1
Dim x As BigInteger = 0
Dim ys As BigInteger = 0
While g = 1
x = y
For i = 1 To r
y = ((y * y) Mod n + c) Mod n
Next
Dim k = New BigInteger(0)
While (k < r And g = 1)
ys = y
For i = 1 To BigInteger.Min(m, r - k)
y = ((y * y) Mod n + c) Mod n
q = q * (BigInteger.Abs(x - y)) Mod n
Next
g = gcd(q, n)
k = k + m
End While
r = r * 2
End While
If g = n Then
While True
ys = ((ys * ys) Mod n + c) Mod n
g = gcd(BigInteger.Abs(x - ys), n)
If g > 1 Then
Exit While
End If
End While
End If
Return g
End Function