求一个数的最大素数的算法

Posted

技术标签:

【中文标题】求一个数的最大素数的算法【英文标题】:Algorithm to find Largest prime factor of a number 【发布时间】:2010-09-06 14:26:12 【问题描述】:

计算一个数的最大质因数的最佳方法是什么?

我认为最有效的方法如下:

    找到能整除的最小素数 检查除法结果是否为素数 如果不是,则查找下一个最低值 转到 2。

我的这个假设是基于更容易计算小的素因数。这是对的吗?我应该研究哪些其他方法?

编辑:我现在意识到,如果有两个以上的质数因子在起作用,我的方法是徒劳的,因为当结果是其他两个质数的乘积时,第 2 步会失败,因此需要递归算法。

再次编辑:现在我意识到这仍然有效,因为最后找到的素数必须是最高的,因此对第 2 步的非素数结果的任何进一步测试都会导致更小的素数.

【问题讨论】:

我的方法是:(1)将大的,可能的数字除以2; (2) 检查大数是否均分; (3) 如果是,检查除以 2 的数是否为素数。如果是,请将其退回。 (4) 否则,将除以 2 的数字减去 1,返回步骤 3。 1. 找到任何可以明确划分的数字(对于 i = 2 到 int(sqr(num)))2. 除以该数字(num = num/i)并重复直到没有发现在1.的区间3.num是最大因子 我们可以用小素数除,最后剩下的就是最大素数(我猜) 【参考方案1】:

我认为最好将所有可能的素数存储在某个地方,然后 n 并遍历它们以找到最大的除数。您可以从prime-numbers.org 获取素数。

当然我假设你的数字不是太大:)

【讨论】:

【参考方案2】:

这可能并不总是更快,但更乐观的是你会找到一个大的素除数:

    N 是你的号码 如果是素数,那么return(N) 计算质数直到Sqrt(N) 按降序遍历素数(最大在前) 如果N is divisible by Prime 那么Return(Prime)

编辑:在第 3 步中,您可以使用 Eratosthenes 筛或阿特金斯筛或任何您喜欢的东西,但筛子本身不会为您找到最大的主要因素。 (这就是为什么我不会选择 SQLMenace 的帖子作为官方答案的原因......)

【讨论】:

您不需要进行试因式分解来确定它是否为质数(步骤 2)吗?另外,考虑找出 15 的最大素数。直到 sqrt(15) 的素数是 2 和 3;但最大的质因数是 5,不是吗?与 20 类似。【参考方案3】:
n = abs(number);
result = 1;
if (n mod 2 == 0) 
 result = 2;
 while (n mod 2 = 0) n /= 2;

for(i=3; i<sqrt(n); i+=2) 
 if (n mod i == 0) 
   result = i;
   while (n mod i = 0)  n /= i;
 

return max(n,result)

有一些模检验是多余的,因为如果所有因素 2 和 3 都已被删除,n 永远不能被 6 整除。您只能允许 i 使用素数,这在此处的其他几个答案中有所显示。

您实际上可以在这里将 Eratosthenes 的筛子交织在一起:

首先创建整数列表 到sqrt(n)。 在 for 循环中标记所有倍数 我直到新的sqrt(n) 素数,并改用 while 循环。 将 i 设置为中的下一个素数 列表。

另见this question。

【讨论】:

【参考方案4】:

实际上有几种更有效的方法可以找到大数的因数(对于较小的数,试除法工作得相当好)。

如果输入数字有两个非常接近其平方根的因子,则一种非常快的方法称为Fermat factorisation。它利用恒等式 N = (a + b)(a - b) = a^2 - b^2,易于理解和实现。不幸的是,它通常不是很快。

最广为人知的分解长达 100 位数字的方法是 Quadratic sieve。作为奖励,部分算法很容易通过并行处理完成。

我听说的另一种算法是Pollard's Rho algorithm。它不如一般的二次筛高效,但似乎更容易实现。


一旦你决定如何将一个数字分成两个因子,这里是我能想到的最快的算法来找到一个数字的最大素因子:

创建一个优先级队列,该队列最初存储数字本身。每次迭代,您从队列中删除最大的数字,并尝试将其分成两个因素(当然,不允许 1 成为这些因素之一)。如果这一步失败,数字是质数,你有你的答案!否则,您将这两个因素添加到队列中并重复。

【讨论】:

Pollard rho 和椭圆曲线法在去除数字的小素因数方面比二次筛法要好得多。无论数量多少,QS 的运行时间都差不多。哪种方法更快取决于您的号码是多少; QS 将更快地破解难以分解的数字,而 rho 和 ECM 将更快地破解易于分解的数字。 感谢您提供二次变化建议。我需要为我的一个客户实现这一点,我提出的初始版本与@mercutio 在他的问题中建议的内容类似。二次解是我现在在math.tools/calculator/prime-factors 为我的客户的工具提供动力的原因。 如果有解决这个问题的有效方法,那是不是意味着网络加密不安全?【参考方案5】:

最简单的解决方案是一对相互递归函数。

第一个函数生成所有素数:

    从所有大于 1 的自然数列表开始。 删除所有非质数。也就是说,没有素因数的数字(除了它们自己)。见下文。

第二个函数以升序返回给定数字n 的质因数。

    列出所有素数(见上文)。 删除所有不是n 的因数的数字。

n 的最大素数是第二个函数给出的最后一个数。

此算法需要一个惰性列表或具有按需调用语义的语言(或数据结构)。

为澄清起见,这里是 Haskell 中上述的一种(低效)实现:

import Control.Monad

-- All the primes
primes = 2 : filter (ap (<=) (head . primeFactors)) [3,5..]

-- Gives the prime factors of its argument
primeFactors = factor primes
  where factor [] n = []
        factor xs@(p:ps) n =
          if p*p > n then [n]
          else let (d,r) = divMod n p in
            if r == 0 then p : factor xs d
            else factor ps n

-- Gives the largest prime factor of its argument
largestFactor = last . primeFactors

加快速度只是更聪明地检测哪些数字是素数和/或n的因数,但算法保持不变。

【讨论】:

【参考方案6】:

所有数都可以表示为素数的乘积,例如:

102 = 2 x 3 x 17
712 = 2 x 2 x 2 x 89

您可以通过简单地从 2 开始并继续除法直到结果不是您的数的倍数来找到这些:

712 / 2 = 356 .. 356 / 2 = 178 .. 178 / 2 = 89 .. 89 / 89 = 1

使用这种方法,您不必实际计算任何素数:它们都将是素数,因为您已经尽可能多地分解了所有前面的数字。

number = 712;
currNum = number;    // the value we'll actually be working with
for (currFactor in 2 .. number) 
    while (currNum % currFactor == 0) 
        // keep on dividing by this number until we can divide no more!
        currNum = currNum / currFactor     // reduce the currNum
    
    if (currNum == 1) return currFactor;    // once it hits 1, we're done.

【讨论】:

是的,但这是非常低效的。一旦你除掉了所有的 2,你真的不应该尝试除以 4,或除以 6,或......;在极限情况下,只检查素数或使用其他算法确实效率更高。 +1 来抵消 wnoise,我认为这是错误的。尝试除以 4 只会发生一次,并且会立即失败。我认为这并不比从候选列表中删除 4 更糟糕,而且肯定比事先找到所有素数要快。 @Beowulf。在投票之前尝试运行此代码。它返回素数;你只是不懂算法。 代码工作正常,但如果传入的数字是素数,则速度很慢。我也只会跑到正方形并增加 2。不过,对于非常大的数字来说,这可能太慢了。 blabla999 完全正确。例如 1234567898766700 = 2*2*5*5*12345678987667。当我们到达currFactor = 3513642 时,我们知道 12345678987667 是素数,并且应该将其作为答案返回。相反,此代码将继续枚举,直到 12345678987667 本身。这比必要的慢 3,513,642 倍。【参考方案7】:

在我看来,给出的算法的第 2 步并不是一种高效的方法。你没有合理的期望它是素数。

此外,先前暗示埃拉托色尼筛法的答案是完全错误的。我刚刚写了两个程序来分解 123456789。一个基于 Sieve,一个基于以下:

1)  Test = 2 
2)  Current = Number to test 
3)  If Current Mod Test = 0 then  
3a)     Current = Current Div Test 
3b)     Largest = Test
3c)     Goto 3. 
4)  Inc(Test) 
5)  If Current < Test goto 4
6)  Return Largest

这个版本比 Sieve 快 90 倍。

问题是,在现代处理器上,操作类型的重要性远低于操作数量,更不用说上面的算法可以在缓存中运行,而 Sieve 则不能。 Sieve 使用了很多操作来剔除所有的合数。

另外请注意,我将确定的因素分开会减少必须测试的空间。

【讨论】:

这就是我所说的,但被否决了:(我想问题是如果这个数字有一个非常大的素因数(比如它自己),那么这个方法必须一直循环到这个数字。但在很多情况下,这种方法非常有效。 读回你的也是一样,但你的第一部分令人困惑。 试试这个号码 143816789988504044536402352738195137863656439,让我知道这有多高效......【参考方案8】:

这是我所知道的最好的算法(在 Python 中)

def prime_factors(n):
    """Returns all the prime factors of a positive integer"""
    factors = []
    d = 2
    while n > 1:
        while n % d == 0:
            factors.append(d)
            n /= d
        d = d + 1

    return factors


pfs = prime_factors(1000)
largest_prime_factor = max(pfs) # The largest element in the prime factor list

上述方法在最坏的情况下(输入为素数时)在O(n) 中运行。

编辑: 下面是O(sqrt(n)) 版本,正如评论中所建议的那样。这是代码,再一次。

def prime_factors(n):
    """Returns all the prime factors of a positive integer"""
    factors = []
    d = 2
    while n > 1:
        while n % d == 0:
            factors.append(d)
            n /= d
        d = d + 1
        if d*d > n:
            if n > 1: factors.append(n)
            break
    return factors


pfs = prime_factors(1000)
largest_prime_factor = max(pfs) # The largest element in the prime factor list

【讨论】:

请在投票之前阅读和/或运行此代码。它工作正常。只需复制和粘贴。正如所写的 prime_factors(1000) 将返回 [2,2,2,5,5,5],它应该被解释为 2^3*5^3,也就是素数分解。 "在最坏的情况下以O(sqrt(n)) 运行" - 不,在最坏的情况下它以O(n) 运行(例如,当n 为素数时。) 容易做到 O(sqrt(n)),当 d*d > n 时停止循环,如果此时 n > 1,那么它的值应该附加到主要因素。 这个有名字吗? 因为 2 是唯一的偶素数,所以不是每次都加 1,您可以单独迭代 d=2,然后将其递增 1,然后从 d=3 开始,您可以递增2.所以它会减少迭代次数... :)【参考方案9】:

我的回答是基于Triptych 的,但在此基础上改进了很多。它基于这样一个事实,即除了 2 和 3,所有质数都是 6n-1 或 6n+1 的形式。

var largestPrimeFactor;
if(n mod 2 == 0)

    largestPrimeFactor = 2;
    n = n / 2 while(n mod 2 == 0);

if(n mod 3 == 0)

    largestPrimeFactor = 3;
    n = n / 3 while(n mod 3 == 0);


multOfSix = 6;
while(multOfSix - 1 <= n)

    if(n mod (multOfSix - 1) == 0)
    
        largestPrimeFactor = multOfSix - 1;
        n = n / largestPrimeFactor while(n mod largestPrimeFactor == 0);
    

    if(n mod (multOfSix + 1) == 0)
    
        largestPrimeFactor = multOfSix + 1;
        n = n / largestPrimeFactor while(n mod largestPrimeFactor == 0);
    
    multOfSix += 6;

我最近写了一个blog article 解释这个算法是如何工作的。

我敢冒险,一种不需要对素性进行测试(也不需要筛子构造)的方法会比使用这些方法运行得更快。如果是这样的话,这可能是这里最快的算法。

【讨论】:

你实际上可以把这个想法更进一步,例如超过 2,3,5 所有素数的形式为 30n+k (n >= 0),其中 k 只取 1 到 29 之间不能被 2,3 或 5 整除的值,即 7,11,13, 17,19,23,29。您甚至可以在迄今为止发现的每几个素数之后动态调整为 2*3*5*7*...*n+k ,其中 k 不能被这些素数中的任何一个整除(请注意,并非所有可能的 k 都需要是素数,例如对于 210n+k,你必须包括 121,否则你会错过331) 我猜应该是while (multOfSix - 1 &lt;= n)【参考方案10】:
#include<stdio.h>
#include<conio.h>
#include<math.h>
#include <time.h>

factor(long int n)

long int i,j;
while(n>=4)
 
if(n%2==0)   n=n/2;   i=2;   

 else
  i=3;
j=0;
  while(j==0)
  
   if(n%i==0)
   j=1;
   n=n/i;
   
   i=i+2;
  
 i-=2;
 
 
return i;
 

 void main()
  
  clock_t start = clock();
  long int n,sp;
  clrscr();
  printf("enter value of n");
  scanf("%ld",&n);
  sp=factor(n);
  printf("largest prime factor is %ld",sp);

  printf("Time elapsed: %f\n", ((double)clock() - start) / CLOCKS_PER_SEC);
  getch();
 

【讨论】:

【参考方案11】:

我知道这不是一个快速的解决方案。发布为希望更容易理解缓慢的解决方案。

 public static long largestPrimeFactor(long n) 

        // largest composite factor must be smaller than sqrt
        long sqrt = (long)Math.ceil(Math.sqrt((double)n));

        long largest = -1;

        for(long i = 2; i <= sqrt; i++) 
            if(n % i == 0) 
                long test = largestPrimeFactor(n/i);
                if(test > largest) 
                    largest = test;
                
            
        

        if(largest != -1) 
            return largest;
        

        // number is prime
        return n;
     

【讨论】:

【参考方案12】:

这里是作为生成器提供的同样的函数@Triptych,它也被稍微简化了。

def primes(n):
    d = 2
    while (n > 1):
        while (n%d==0):
            yield d
            n /= d
        d += 1

然后可以使用以下方法找到最大素数:

n= 373764623
max(primes(n))

以及使用以下方法找到的因素列表:

list(primes(n))

【讨论】:

【参考方案13】:

首先计算一个存储素数的列表,例如2 3 5 7 11 13 ...

每次对一个数进行素数分解时,请使用 Triptych 的实现,但要迭代这个素数列表而不是自然整数。

【讨论】:

【参考方案14】:

这是我在 C# 中的尝试。最后的打印输出是数字的最大素数。我检查了,它可以工作。

namespace Problem_Prime

  class Program
  
    static void Main(string[] args)
    
      /*
       The prime factors of 13195 are 5, 7, 13 and 29.

      What is the largest prime factor of the number 600851475143 ?
       */
      long x = 600851475143;
      long y = 2;
      while (y < x)
      
        if (x % y == 0)
        
          // y is a factor of x, but is it prime
          if (IsPrime(y))
          
            Console.WriteLine(y);
          
          x /= y;
        

        y++;

      
      Console.WriteLine(y);
      Console.ReadLine();
    
    static bool IsPrime(long number)
    
      //check for evenness
      if (number % 2 == 0)
      
        if (number == 2)
        
          return true;
        
        return false;
      
      //don't need to check past the square root
      long max = (long)Math.Sqrt(number);
      for (int i = 3; i <= max; i += 2)
      
        if ((number % i) == 0)
        
          return false;
        
      
      return true;
    

  

【讨论】:

【参考方案15】:

使用 Java:

对于int 值:

public static int[] primeFactors(int value) 
    int[] a = new int[31];
    int i = 0, j;
    int num = value;
    while (num % 2 == 0) 
        a[i++] = 2;
        num /= 2;
    
    j = 3;
    while (j <= Math.sqrt(num) + 1) 
        if (num % j == 0) 
            a[i++] = j;
            num /= j;
         else 
            j += 2;
        
    
    if (num > 1) 
        a[i++] = num;
    
    int[] b = Arrays.copyOf(a, i);
    return b;

对于long 值:

static long[] getFactors(long value) 
    long[] a = new long[63];
    int i = 0;
    long num = value;
    while (num % 2 == 0) 
        a[i++] = 2;
        num /= 2;
    
    long j = 3;
    while (j <= Math.sqrt(num) + 1) 
        if (num % j == 0) 
            a[i++] = j;
            num /= j;
         else 
            j += 2;
        
    
    if (num > 1) 
        a[i++] = num;
    
    long[] b = Arrays.copyOf(a, i);
    return b;

【讨论】:

【参考方案16】:
    //this method skips unnecessary trial divisions and makes 
    //trial division more feasible for finding large primes

    public static void main(String[] args) 
    
        long n= 1000000000039L; //this is a large prime number 
        long i = 2L;
        int test = 0;

        while (n > 1)
        
            while (n % i == 0)
            
                n /= i;     
            

            i++;

            if(i*i > n && n > 1) 
            
                System.out.println(n); //prints n if it's prime
                test = 1;
                break;
            
        

        if (test == 0)  
            System.out.println(i-1); //prints n if it's the largest prime factor
    

【讨论】:

您是否尝试过使用 1,000,000,000,039 的代码?它也应该在眨眼间运行。是吗? 你可以提前知道,无需尝试。 10^12 = (2*5)^12 = 2^12 * 5^12。所以你的while 循环将通过i2,2,2,2,2,2,2,2,2,2,2,2, 2,3,4,5, 2,3,4,5, 2,3,4,5, 2,3,4,5, 2,3,4,5, 2,3,4,5, 2,3,4,5, 2,3,4,5, 2,3,4,5, 2,3,4,5, 2,3,4,5, 2,3,4,5 值。全部 60 次迭代。但是对于 (10^12+39) 将有 (10^12+38) 次迭代,i=2,3,4,5,6,...,10^12+39。即使 10^10 次操作需要一秒钟,10^12 也需要 100 秒。但实际上只需要 10^6 次迭代,如果 10^10 次操作需要 1 秒,10^6 将需要 1/10,000 秒。 因为我没有意识到 (10^12+39) 是我现在所做的素数。我明白你在说什么。 好的,所以你可以改变你的代码,这样素数就不会有这么大的减速:如果 n = ab 和 a a a = n,即 aa long n = 2*1000000000039L 时会发生什么?它是否以应有的速度工作? (另外,你能通过使用return; 语句来简化你的代码吗?)。 (如果你想让我停止轻推你,就这么说吧;))【参考方案17】:
#python implementation
import math
n = 600851475143
i = 2
factors=set([])
while i<math.sqrt(n):
   while n%i==0:
       n=n/i
       factors.add(i)
   i+=1
factors.add(n)
largest=max(factors)
print factors
print largest

【讨论】:

25 是 25 的最大素数吗?【参考方案18】:

在 C++ 中使用递归计算数字的最大素数。代码的工作原理解释如下:

int getLargestPrime(int number) 
    int factor = number; // assumes that the largest prime factor is the number itself
    for (int i = 2; (i*i) <= number; i++)  // iterates to the square root of the number till it finds the first(smallest) factor
        if (number % i == 0)  // checks if the current number(i) is a factor
            factor = max(i, number / i); // stores the larger number among the factors
            break; // breaks the loop on when a factor is found
        
    
    if (factor == number) // base case of recursion
        return number;
    return getLargestPrime(factor); // recursively calls itself

【讨论】:

【参考方案19】:

通过从数字中删除所有质因数的 Python 迭代方法

def primef(n):
    if n <= 3:
        return n
    if n % 2 == 0:
        return primef(n/2)
    elif n % 3 ==0:
        return primef(n/3)
    else:
        for i in range(5, int((n)**0.5) + 1, 6):
            #print i
            if n % i == 0:
                return primef(n/i)
            if n % (i + 2) == 0:
                return primef(n/(i+2))
    return n

【讨论】:

【参考方案20】:

这是我快速计算最大素数的方法。 它基于修改后的x 不包含非主要因素的事实。为了实现这一点,我们会在找到一个因子后立即划分x。然后,唯一剩下的就是返回最大的因子。应该已经是素数了。

代码(Haskell):

f max' x i | i > x = max'
           | x `rem` i == 0 = f i (x `div` i) i  -- Divide x by its factor
           | otherwise = f max' x (i + 1)  -- Check for the next possible factor

g x = f 2 x 2

【讨论】:

但这不会也尝试除以所有偶数吗?【参考方案21】:

javascript 代码:

'option strict';

function largestPrimeFactor(val, divisor = 2)  
    let square = (val) => Math.pow(val, 2);

    while ((val % divisor) != 0 && square(divisor) <= val) 
        divisor++;
    

    return square(divisor) <= val
        ? largestPrimeFactor(val / divisor, divisor)
        : val;

用法示例:

let result = largestPrimeFactor(600851475143);

Here is an example of the code:

【讨论】:

【参考方案22】:

以下 C++ 算法不是最好的,但它适用于十亿以下的数字并且速度非常快

#include <iostream>
using namespace std;

// ------ is_prime ------
// Determines if the integer accepted is prime or not
bool is_prime(int n)
    int i,count=0;
    if(n==1 || n==2)
      return true;
    if(n%2==0)
      return false;
    for(i=1;i<=n;i++)
    if(n%i==0)
        count++;
    
    if(count==2)
      return true;
    else
      return false;
 
 // ------ nextPrime -------
 // Finds and returns the next prime number
 int nextPrime(int prime)
     bool a = false;
     while (a == false)
         prime++;
         if (is_prime(prime))
            a = true;
     
  return prime;
 
 // ----- M A I N ------
 int main()

      int value = 13195;
      int prime = 2;
      bool done = false;

      while (done == false)
          if (value%prime == 0)
             value = value/prime;
             if (is_prime(value))
                 done = true;
             
           else 
             prime = nextPrime(prime);
          
      
        cout << "Largest prime factor: " << value << endl;
 

【讨论】:

【参考方案23】:

我使用的算法继续将数字除以当前的素数。

我在 python 3 中的解决方案:

def PrimeFactor(n):
    m = n
    while n%2==0:
        n = n//2
    if n == 1:         # check if only 2 is largest Prime Factor 
        return 2
    i = 3
    sqrt = int(m**(0.5))  # loop till square root of number
    last = 0              # to store last prime Factor i.e. Largest Prime Factor
    while i <= sqrt :
        while n%i == 0:
            n = n//i       # reduce the number by dividing it by it's Prime Factor
            last = i
        i+=2
    if n> last:            # the remaining number(n) is also Factor of number 
        return n
    else:
        return last
print(PrimeFactor(int(input()))) 

输入:10 输出:5

输入:600851475143 输出:6857

【讨论】:

【参考方案24】:

类似于@Triptych 的答案,但也不同。在此示例中,未使用列表或字典。代码是用 Ruby 编写的

def largest_prime_factor(number)
  i = 2
  while number > 1
    if number % i == 0
      number /= i;
    else
      i += 1
    end
  end
  return i
end

largest_prime_factor(600851475143)
# => 6857

【讨论】:

终于有了可读性和即时(在 js 中)同时执行的东西。我正在尝试使用无限素数列表,但它在 100 万时已经太慢了。【参考方案25】:

“James Wang”在网上找到了这个解决方案

public static int getLargestPrime( int number) 

    if (number <= 1) return -1;

    for (int i = number - 1; i > 1; i--) 
        if (number % i == 0) 
            number = i;
        
    
    return number;

【讨论】:

【参考方案26】:

使用筛子的素因子:

#include <bits/stdc++.h>
using namespace std;
#define N 10001  
typedef long long ll;
bool visit[N];
vector<int> prime;

void sieve()

            memset( visit , 0 , sizeof(visit));
            for( int i=2;i<N;i++ )
            
                if( visit[i] == 0)
                
                    prime.push_back(i);
                    for( int j=i*2; j<N; j=j+i )
                    
                        visit[j] = 1;
                    
                
               

void sol(long long n, vector<int>&prime)

            ll ans = n;
            for(int i=0; i<prime.size() || prime[i]>n; i++)
            
                while(n%prime[i]==0)
                
                    n=n/prime[i];
                    ans = prime[i];
                
            
            ans = max(ans, n);
            cout<<ans<<endl;

int main() 

           ll tc, n;
           sieve();

           cin>>n;
           sol(n, prime);

           return 0;

【讨论】:

【参考方案27】:

这是我在 Clojure 中的尝试。只为prime? 和素数的素数走赔率,即。 sieve。使用惰性序列有助于在需要它们之前生成值。

(defn prime? 
  ([n]
    (let [oddNums (iterate #(+ % 2) 3)]
    (prime? n (cons 2 oddNums))))
  ([n [i & is]]
    (let [q (quot n i)
          r (mod n i)]
    (cond (< n 2)       false
          (zero? r)     false
          (> (* i i) n) true
          :else         (recur n is)))))

(def primes 
  (let [oddNums (iterate #(+ % 2) 3)]
  (lazy-seq (cons 2 (filter prime? oddNums)))))

;; Sieve of Eratosthenes
(defn sieve
  ([n] 
    (sieve primes n))
  ([[i & is :as ps] n]
    (let [q (quot n i)
          r (mod n i)]
    (cond (< n 2)       nil
          (zero? r)     (lazy-seq (cons i (sieve ps q)))
          (> (* i i) n) (when (> n 1) (lazy-seq [n]))
          :else         (recur is n)))))

(defn max-prime-factor [n]
  (last (sieve n)))

【讨论】:

【参考方案28】:

猜猜,除了执行分解之外没有直接的方法,就像上面的例子所做的那样,即

在迭代中,您确定数字 N 的“小”因子 f,然后继续简化问题“找到 N 的最大素因子” :=N/f 与候选因子 >=f ".

f 的特定大小开始,如果您对减少的 N' 进行素性测试,则预期的搜索时间会更短,以防万一,您的 N'已经是初始N的最大素数了。

【讨论】:

【参考方案29】:

受您的问题启发,我决定在 Python 中实现我自己的因式分解(并找到最大的素因数)。

据我所知,可能最容易实现但非常有效的因式分解算法是Pollard's Rho 算法。它的运行时间最多为O(N^(1/4)),比试除算法的时间O(N^(1/2))快得多。两种算法只有在复合(非素数)数的情况下才有这些运行时间,这就是为什么应该使用素数测试来过滤掉素数(不可因数)。

我在代码中使用了以下算法:Fermat Primality Test ...、Pollard's Rho Algorithm ...、Trial Division Algorithm。在运行 Pollard 的 Rho 之前使用 Fermat 素数测试以过滤掉素数。 Trial Division 被用作后备,因为 Pollard 的 Rho 在极少数情况下可能无法找到因子,尤其是对于一些小数字。

显然,在将一个数字完全分解为已排序的素因子列表后,最大的素因子将是该列表中的最后一个元素。在一般情况下(对于任何随机数),除了完全分解一个数字之外,我不知道有任何其他方法可以找出最大的素因数。

作为我的代码中的示例,我将 Pi 的第一个 190 小数位分解,代码在 1 秒内分解这个数字,并显示最大的素数因子,即 165位数(545 位)!

Try it online!

def is_fermat_probable_prime(n, *, trials = 32):
    # https://en.wikipedia.org/wiki/Fermat_primality_test
    import random
    if n <= 16:
        return n in (2, 3, 5, 7, 11, 13)
    for i in range(trials):
        if pow(random.randint(2, n - 2), n - 1, n) != 1:
            return False
    return True

def pollard_rho_factor(N, *, trials = 16):
    # https://en.wikipedia.org/wiki/Pollard%27s_rho_algorithm
    import random, math
    for j in range(trials):
        i, stage, y, x = 0, 2, 1, random.randint(1, N - 2)
        while True:
            r = math.gcd(N, x - y)
            if r != 1:
                break
            if i == stage:
                y = x
                stage <<= 1
            x = (x * x + 1) % N
            i += 1
        if r != N:
            return [r, N // r]
    return [N] # Pollard-Rho failed

def trial_division_factor(n, *, limit = None):
    # https://en.wikipedia.org/wiki/Trial_division
    fs = []
    while n & 1 == 0:
        fs.append(2)
        n >>= 1
    d = 3
    while d * d <= n and limit is None or d <= limit:
        q, r = divmod(n, d)
        if r == 0:
            fs.append(d)
            n = q
        else:
            d += 2
    if n > 1:
        fs.append(n)
    return fs

def factor(n):
    if n <= 1:
        return []
    if is_fermat_probable_prime(n):
        return [n]
    fs = trial_division_factor(n, limit = 1 << 12)
    if len(fs) >= 2:
        return sorted(fs[:-1] + factor(fs[-1]))
    fs = pollard_rho_factor(n)
    if len(fs) >= 2:
        return sorted([e1 for e0 in fs for e1 in factor(e0)])
    return trial_division_factor(n)

def demo():
    import time, math
    # http://www.math.com/tables/constants/pi.htm
    # pi = 3.
    #     1415926535 8979323846 2643383279 5028841971 6939937510 5820974944 5923078164 0628620899 8628034825 3421170679
    #     8214808651 3282306647 0938446095 5058223172 5359408128 4811174502 8410270193 8521105559 6446229489 5493038196
    # n = first 190 fractional digits of Pi
    n =   1415926535_8979323846_2643383279_5028841971_6939937510_5820974944_5923078164_0628620899_8628034825_3421170679_8214808651_3282306647_0938446095_5058223172_5359408128_4811174502_8410270193_8521105559_6446229489
    print('Number:', n)
    tb = time.time()
    fs = factor(n)
    print('All Prime Factors:', fs)
    print('Largest Prime Factor:', f'(math.log2(fs[-1]):.02f bits, len(str(fs[-1])) digits)', fs[-1])
    print('Time Elapsed:', round(time.time() - tb, 3), 'sec')

if __name__ == '__main__':
    demo()

输出:

Number: 1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679821480865132823066470938446095505822317253594081284811174502841027019385211055596446229489
All Prime Factors: [3, 71, 1063541, 153422959, 332958319, 122356390229851897378935483485536580757336676443481705501726535578690975860555141829117483263572548187951860901335596150415443615382488933330968669408906073630300473]
Largest Prime Factor: (545.09 bits, 165 digits) 122356390229851897378935483485536580757336676443481705501726535578690975860555141829117483263572548187951860901335596150415443615382488933330968669408906073630300473
Time Elapsed: 0.593 sec

【讨论】:

以上是关于求一个数的最大素数的算法的主要内容,如果未能解决你的问题,请参考以下文章

AtCoder Beginner Contest 142D题判断素数的模板+求一个数的因子的模板

素数定理-欧几里得算法-乘法逆元

算法学习——递推之超级素数

C语言求范围之内素数的和与个数?

LightOJ1356 最大独立集 HK算法 素数分解

两道经典算法题-吉比特2017秋招笔试