大数的有效素数分解

Question

我一直在解决一个小问题，我需要将 18 位数字计算到它们各自的素数分解中。考虑到它确实有效，一切都可以编译并且运行得很好，但我希望减少素数分解的运行时间。我已经实现了递归和线程，但我认为我可能需要一些帮助来理解可能的大量计算算法。

每次我对预先制作的 4 个数字运行此操作时，大约需要 10 秒。如果有任何想法，我想将其减少到可能的 0.06 秒。

我注意到一些算法，例如 Sieve of Eratosthenes，并在计算之前生成了所有素数的列表。我只是想知道是否有人可以详细说明。例如，我在理解如何在我的程序中实施埃拉托色尼筛法或者它是否是一个好主意时遇到了问题。关于如何更好地解决这个问题的任何和所有指示都会非常有帮助！

这是我的代码：

#include <iostream>
#include <thread>
#include <vector>
#include <chrono>

using namespace std;
using namespace std::chrono;

vector<thread> threads;
vector<long long> inputVector;
bool developer = false; 
vector<unsigned long long> factor_base;
vector<long long> primeVector;

class PrimeNumber
{
    long long initValue;        // the number being prime factored
    vector<long long> factors;  // all of the factor values
public:
    void setInitValue(long long n)
    {
        initValue = n;
    }
    void addToVector(long long m)
    {
        factors.push_back(m);
    }
    void setVector(vector<long long> m)
    {
        factors = m;
    }
    long long getInitValue()
    {
        return initValue;
    }
    vector<long long> getVector()
    {
        return factors;
    }
};

vector<PrimeNumber> primes;

// find primes recursively and have them returned in vectors
vector<long long> getPrimes(long long n, vector<long long> vec)
{
    double sqrt_of_n = sqrt(n);

    for (int i = 2; i <= sqrt_of_n; i++)
    {
        if (n % i == 0) 
        {
            return vec.push_back(i), getPrimes(n / i, vec); //cause recursion
        }
    }

    // pick up the last prime factorization number
    vec.push_back(n);

    //return the finished vector
    return vec;
}

void getUserInput()
{
    long long input = -1;
    cout << "Enter all of the numbers to find their prime factors. Enter 0 to compute" << endl;
    do
    {
        cin >> input;
        if (input == 0)
        {
            break;
        }
        inputVector.push_back(input);
    } while (input != 0);
}

int main() 
{

    vector<long long> temp1;   // empty vector
    vector<long long> result1; // temp vector

    if (developer == false)
    {
        getUserInput();
    }
    else
    {
        cout << "developer mode active" << endl;
        long long a1 = 771895004973090566;
        long long b1 = 788380500764597944;
        long long a2 = 100020000004324000;
        long long b2 = 200023423420000000;
        inputVector.push_back(a1);
        inputVector.push_back(b2);
        inputVector.push_back(b1);
        inputVector.push_back(a2);
    }

    high_resolution_clock::time_point time1 = high_resolution_clock::now();

    // give each thread a number to comput within the recursive function
    for (int i = 0; i < inputVector.size(); i++)
    {   
        PrimeNumber prime;
        prime.setInitValue(inputVector.at(i));
        threads.push_back(thread([&]{
            prime.setVector(result1 = getPrimes(inputVector.at(i), temp1));
            primes.push_back(prime);
        }));
    }

    // allow all of the threads to join back together.
    for (auto& th : threads)
    {
        cout << th.get_id() << endl;
        th.join();
    }

    high_resolution_clock::time_point time2 = high_resolution_clock::now();

    // print all of the information
    for (int i = 0; i < primes.size(); i++)
    {
        vector<long long> temp = primes.at(i).getVector();

        for (int m = 0; m < temp.size(); m++)
        {
            cout << temp.at(m) << " ";
        }
        cout << endl;
    }

    cout << endl;

    // so the running time
    auto duration = duration_cast<microseconds>(time2 - time1).count();

    cout << "Duration: " << (duration / 1000000.0) << endl;

    return 0;
}

Answer 1

试除法只适用于分解小数。对于高达 2^64 的 n，您将需要一个更好的算法：我建议从轮分解开始以获取小因子，然后使用 Pollard 的 rho 算法来获取其余部分。试除法是 O(sqrt(n))，rho 是 O(sqrt(sqrt(n)))，所以要快得多。对于 2^64，sqrt(n) = 2^32，但是 sqrt(sqrt(n)) = 2^16，这是一个巨大的改进。您应该期望最多在几毫秒内计算出您的数字。

我没有用于分解的 C++ 代码，但我有可读的 Python 代码。如果你想让我发布它，请告诉我。如果您想了解更多关于车轮分解和 rho 算法的信息，我在my blog 有很多素数资料。

Answer 2

我发现现代处理器上的埃拉托色尼筛法会破坏缓存，因此主内存带宽是限制因素。我在尝试运行多个线程并且未能达到我希望的速度时发现了这一点。

因此，我建议将筛子分成适合 L3 缓存的段。此外，如果从位向量中排除 2、3 和 5 的倍数，则 8 位字节可以表示数轴上的 30 个数字，每个数字 1、7、11、13、17、19 , 23 或 29 以 30 为模——因此，10^9 以内的素数的位图需要 ~32MB——10^9 / (30 * 1024 * 1024)。这几乎是仅排除 2 的倍数的位图大小的一半，即 ~60MB -- 10^9 / (2 * 8 * 1024 * 1024)。

显然，要将筛子运行到 10^9，您需要素数到 sqrt(10^9) - 这需要大约 1,055 个字节，您可以从中生成完整筛子的任何部分，最多 10^9 .

FWIW，我在适中的 AMD Phenom II x6 1090T（8MB L3 缓存）上得到的结果，对于最高 10^9 的素数是：

  1. 1 core,   1 segment    3.260 seconds elapsed
  2. 5 cores,  1 segment    1.830 seconds elapsed
  3. 1 core,   8 segments   1.800 seconds elapsed
  4. 5 cores, 40 segments   0.370 seconds elapsed

我所说的“段”是指筛子的一部分。在这种情况下，筛子约为 32MB，因此在有多个段的情况下，它们在任何时候都使用大约 4MB 的 L3 缓存。

这些时间包括扫描完成的筛子并将所有素数生成为整数数组所需的时间。这需要大约 0.5 秒的 CPU 时间！因此，在不实际从中提取素数的情况下运行筛子，在上述情况 (4) 中需要 0.270 秒。

FWIW，通过使用预先计算的模式初始化每个段，删除 7、11、13 和 17 的倍数，我得到了一个小的改进 - 在情况 (4) 中为 0.240 秒。该模式是 17,017 字节。

显然，要在 0.06 秒内进行一次分解...您需要预先计算筛子！

Answer 3

for(int i  = 2; i * i <= n; ++i) //no sqrt, please
{
    while(n%i == 0) //while, not if
    {
         factors.push_back(i);
         n/=i;
    }
}
if(n != 1)
{
    factors.push_back(n);
}

这基本上是您算法的一个更简洁的实现。它的复杂度是 N 的平方。即使是 18 位数字，它也能很快工作，但前提是主要因素都很小。如果它是两个大素数的乘积，或者更糟糕的是，它本身就是素数，这将运行大约 10 秒。

Answer 4

通过更改循环可以轻松实现 2 的简单加速：

if (n % 2) {
    return vec.push_back(i), getPrimes(n / i, vec);
}

for (int i = 3; i <= sqrt_of_n; i += 2)
{
    if (n % i == 0) 
    {
        return vec.push_back(i), getPrimes(n / i, vec); //cause recursion
    }
}

你首先应该用 2 来测试这个数字。然后，从 3 开始，您再次测试，一次将循环增加两个。您已经知道 4, 6, 8, ... 是偶数并且有 2 作为因数。对偶数进行测试可以将复杂性降低一半。

要分解一个数字N，您只需要素数 1e9 的所有素数进行测试，并且由于 there are 98 millon primes 小于 2e9，您可以轻松地在当今的计算机上存储 100 万个数字并并行运行因式分解。如果每个数字占用 8 个字节的 RAM (int64_t)，则 100 个质数将占用 800 MB 的内存。这个算法是SPOJ problem #2, Prime Generator的经典解法。

列出所有适合 32 位 int 的小素数的最佳方法是构建 Eratostenes 筛。我告诉过你，我们需要小于 sqrt(N) 的素数来分解任何 N，所以要分解 64 位整数，你需要所有适合 32 位数字的素数。

Answer 5

算法：

1. 将数字除以 2 直到它不能被它整除，然后存储结果并 显示它。
2. 将数字除以 3 直到它不能被 3 整除并显示结果，
3. 对 5,7... 等重复相同的过程，直到 n 的平方根。

4. 如果最终结果数是素数，则显示其计数为 1。

   int main() {
   long long n;
   cin >> n;
   int count =0 ;
   while(!(n%2)){
       n = n / 2;
       count++;
   }
   if(count > 0) {
       cout<<"2^"<<count<<" ";
   }
   for(long long i=3;i<=sqrt(n) ; i+=2){
       count=0;
       while(n%i == 0){
           count++;
           n = n/i;
       }
       if(count){
           cout << i <<"^" <<count<<" ";
       }
   
   }
   if(n>2){
       cout<<n <<"^1";
   }
   
   
   
   return 0;
   }
   

输入：100000000 输出 2^8 5^8

Answer 6

整数分解算法，非常简单，甚至可以在算盘上实现。

void start()
{ 
   int a=4252361;    //  integer to factorize
   int b=1,  c,  d=a-1;  

   while ((a > b) && (d != b))
  {
    if (d > b)
     {
       c=c+b;
       a=a-1;
       d=a-c-1;
     }
    if (d < b)
     {
       c=b-d;
       b=b+1;
       a=a-1;
       d=a-c-1;
     }         
  }
    if ((d == b)&&(a > b)) 
  Alert ("a = ",  a-1,  " * ", b+1); 

    if ((d < b)&&(a <= b)) 
  Alert ("a  is a prime");

  return;
}

该算法由我 Miliaev Viktor Konstantinovich 编写，他出生于 1950 年 7 月 26 日。 它是用 MQL4 15.12 编写的。 2017. 邮箱：tenfacet27@gmail.com