更好的算法-下一个半素

Question

给定 n，求 m 使得 m 是大于 n 的最小半素数。

下一个素数相当简单，而半素数则不那么简单。需要明确的是，只需要半素数，尽管同时获得因子会很方便。

我想到了一些方法，但我相信有更好的方法。

假设算术运算为O（1）。我使用了埃拉托斯特尼筛，它是O（n log log n），我知道阿特金筛，但我喜欢我的半优化埃拉托斯特尼。

从n开始计数，当你找到一个半素数时停止。

这看起来很愚蠢，但如果有一个O（logn）半素数测试或O（1）测试给定以下素数，这可能比其他2更快。

半素数分布似乎比素数分布高得多，因此通过良好的半素数检验，这实际上可能比 O（n） 更好。

定义prev(x)和next(x ),分别给出前一个和下一个素数，如果素数存储在树中或用列表二分搜索法存储，则可以是O(log n)。

先做筛子。

从 p= prev（平方（n））和 q= 下一个（n/p）开始。而 pq

这可以保证找到正确的答案，但速度很慢。但仍然是O（n log n），因此可能可以接受。

像往常一样从筛子开始。为 O（1） 素数检验创建筛子的哈希集视图。

从 p=2 开始。迭代素数直至平方（n）。对于每个 p，获取 q=（（（无/页 1）/2）*2） 1=（（（无/页 1）

我在Java中实现了#1和#3，两者都使用相同的Eratosthenes实现筛子。大部分运行时间都花在筛选上，所以如果有优化要做，它就在筛子里。经过一些优化后，计数（#1）比质数计数（#3）快一倍，在最后一个也是最大的测试中（11位十进制数字n）。

然而，仍然有希望，因为筛子需要延伸多远取决于质数测试的最大数字。如果存在具有较低素数测试边界的半素数测试，则计数方法可能被证明是更快的。

肯定有更好的算法吗？或者至少有一种更好的方法来实现这一点？

Answer 1

根据@DanielWagner的建议提出的评论（现已删除），这里有一个非优化的半螺旋筛，每个条目使用两个位来保持因子计数。

筛选条目包含发现的因子数量，限于3个。一次标记通过会使相关筛目达到饱和增量。

因为我们也关心两个相等的因子，所以我们也筛分素数的平方。素数的幂可以在筛分过程中识别，因为它们的因子计数将为1（素数计数为0，半素数为2，其他整数为3）。当我们标记素数的平方（这将是遇到的素数的一次幂）时，我们可以对每个条目进行2的饱和加法，但作为微优化，代码只是将计数直接设置为3。

假设筛子不包含偶数的条目（像往常一样），我们将半素数4和所有因子为2和奇数素数的半素数特殊情况。

下面的代码用（伪）C表示，只显示如何进行筛分操作。一些细节，包括saturating_increment和其他筛选访问函数的定义，都被省略了，因为它们很明显，而且只会分散注意力。

/* Call `handler` with every semiprime less than `n`, in order */ 
void semi_sieve(uint32_t n, void(*handler)(uint32_t semi)) {
  Sieve sieve(n);
  if (n > 4) handler(4); /* Special case */
  for (uint32_p p = 3; p < n; p += 2) {
    switch (sieve.get(p)) {
      case 0: /* A prime */
        for (uint32_t j = p + p + p; j < n; j += p + p)
          sieve.saturating_increment(j);
        break;
      case 1: /* The square of a prime */
        handler(p);
        for (uint32_t j = p + p + p; j < n; j += p + p)
          sieve.set(j, 3);
        break;
      case 2: /* A semiprime */
        handler(p);
        break;
      case 3: /* Composite non-semiprime */
        break;
      default: /* Logic error */
    }
    /* If the next number might be twice an odd prime, check the sieve */
    if (p + 1 < n && p % 4 == 1 && 0 == sieve.get((p + 1)/2))
      handler(p + 1);
  }
}

注意：我知道上面的扫描使用了整个范围的质数，而不仅仅是平方根。这必须付出一定的代价，但我认为这只是常数的变化。有可能提前终止扫描，恢复一点常数。

Answer 2

以下是基于我上面的评论的一些代码：我们将运行埃拉托斯特尼的筛子，但在我们这样做时，存储一些额外的数据，而不仅仅是“质数与否”。它是在Haskell中，我认识到它不是最常见的语言，所以我将内联评论每个位的作用。首先，一些库导入：

import Control.Monad
import Control.Monad.ST
import Data.Array.ST
import Data.Maybe

我们将定义一个新的类型，< code>Primality，我们将使用它来存储每个数字的两个质因数。

data Primality
    = Prime
    | OneFactor Integer
    | TwoFactor Integer Integer
    | ManyFactor
    deriving (Eq, Ord, Read, Show)

这表示有四种类型的< code >素性值:或者是值< code >素，对于某个无界整数< code>n是形式为< code>OneFactor n的值，对于两个无界整数< code>n和< code>n'是形式为< code > two factor n ' 的值，或者是值< code>ManyFactor。所以这有点像一个< code>Integer的列表，它最多有两个整数长(或者一个注释说它有三个整数长或更长)。我们可以像这样在列表中添加一些因素:

addFactor :: Integer -> Primality -> Primality
addFactor f Prime = OneFactor f
addFactor f (OneFactor f') = TwoFactor f f'
addFactor _ _ = ManyFactor

给定一个数的素数因子列表，很容易检查它是否是半素数：它最多必须有两个较小的素数因子，其乘积就是数字本身。

isSemiprime :: Integer -> Primality -> Bool
isSemiprime n (OneFactor f   ) = f * f  == n
isSemiprime n (TwoFactor f f') = f * f' == n
isSemiprime _ _ = False

现在我们来写筛子。根据伯特兰的假设，对于任何n，在n/2和n之间有一个素数;这意味着在n和2n之间有一个半素数（即假设给我们的素数是两倍）。更重要的是，任何这样的半素数都不能有一个大于n的因子（从那时起，另一个因子必须小于2！因此，我们将用高达 n 的因子筛选最大 2n 的数字，然后检查 n 和 2n 之间的数字是否存在半素数。因为后一个检查是O（1），所以我们属于你提出的第一种情况。所以：

nextSemiprime :: Integer -> Integer
nextSemiprime n = runST $ do

在2和2n之间创建一个索引数组，在每个位置初始化为Prime。

    arr <- newSTArray (2,2*n) Prime

对于介于 2 和 n 之间的每个潜在素数 p...

    forM_ [2..n] $ \p -> do

…我们目前认为是Prime…

        primality <- readArray arr p
        when (primality == Prime) $

…将p添加到每个p倍数的因子列表中。

            forM_ [2*p,3*p..2*n] $ \i ->
                modifyArray arr i (addFactor p)

现在对n 1和2n之间的剩余数字进行半素数的线性搜索。每次调用isSemiPrime都会花费一次乘法，因此它们是O（1）。从技术上讲，搜索可能会失败；fromJust的

    fromJust <$> findM (\i -> isSemiprime i <$> readArray arr i) [n+1..2*n]

这就是下一个半总理的全部身体。它使用了一些帮助器函数，这些函数确实应该在标准库的某个地方。第一种是线性搜索算法;它只是走下一个列表，寻找一个满足谓词的元素。

findM :: Monad m => (a -> m Bool) -> [a] -> m (Maybe a)
findM f [] = return Nothing
findM f (x:xs) = do
    done <- f x
    if done then return (Just x) else findM f xs

修改Array函数只读取一个数组元素并写回修改后的版本。想想 Arr[ix] = f（arr[ix]）; 在 C 中。

modifyArray :: (MArray a e m, Ix i) => a i e -> i -> (e -> e) -> m ()
modifyArray arr ix f = readArray arr ix >>= writeArray arr ix . f

而newSTArray是必需的，因为Haskell对数组的处理方式变化无常：所有数组操作都是针对您使用的数组类型的多态操作，这既方便又烦人。这告诉编译器我们想要这个程序使用哪种数组。

newSTArray :: Ix i => (i,i) -> e -> ST s (STArray s i e)
newSTArray = newArray

您可以在这里尝试一下，其中包括一个简单的main，用于打印前100个半素数。（尽管后一项任务可以以更有效的方式完成，如果这是目标的话！）

尽管当前的算法只返回下一个半素，但很容易修改它以返回下一个半素的因式分解:只返回相关的< code >素性值，而不是< code >整数本身。

Answer 3

人们回答的问题略有不同，所以让我把它分成几个部分。

1. 是半质数（n）

给定一个值n，它是半素数吗？对于微小的输入，我们当然可以预先计算并返回O（1）或搜索的答案。在某个时候，我们会被存储需求淹没。据我所知，这个问题没有非常有效的解决方案。它类似于素数或无平方测试，因为我们可以通过简单的整除测试快速剔除大多数随机输入。假设我们有一个快速素数测试，其中包括一些预测试，大部分工作只是寻找一个小因子，然后返回余数是否是素数。对于没有小因子的数字，我们可以做一些因数分解（例如：Brent/Pollard Rho）或试验除以n^（1/3）。

在我的Macbook上，对于1e8到1e7 1e7的范围，每个数字大约需要0.4微秒，对于1e16到1e16 1e7的范围，每个数字大约需要2微秒。

对于大的半质数或近似半质数，我不确定有比找到一个单一因子更好的解决方案。我们要求试除法只适用于N^(1/3)但是有更有效的标准因式分解算法。一些实现包括Charles Greathouse、mine和RosettaCode的许多实现。

在1e16时，到下一个半素数的平均距离低于10，很少超过100。像以前一样，如果您想进行预计算，使用内存，并且可以忽略或摊销设置时间，则可以快速回答。但是，一旦过去的小输入，这又变得非常麻烦。

我不认为你能在一个简单的< code > while(1){ n；if (is_semiprime(n))返回n；}假设一个好的is_semiprime例程。做一个完整的筛选对我来说要慢得多，但你的里程数可能会有所不同。一旦输入超过25位数，它就真的无法运行了。您可以通过使用质数幂增加因子计数来进行部分筛选，这意味着我们只需要对明显不是半质数的结果进行完整测试。对我来说没有节省多少时间，这是有意义的，因为我们只删除了一些本机modulos。如果我们看到的是1000位数的输入，那么我认为部分筛选很有意义。

在我的Macbook上，从1e8开始，使用简单的is_semiprime方法连续1e6次调用next_semiprime，每次调用大约需要2微秒，而从1e16开始，每次调用需要17微秒。

一些答案似乎正在考虑这个问题。特别是当低

注意：一个写得好的SoE比SoA快，所以不要被推荐《阿特金筛》的人分心，因为他们可能刚刚阅读了维基百科页面的前几段。当然，筛选、素性测试和预测试的实施细节将对结论产生影响。以及缓存数据的预期输入大小、模式和容差。