查找最小N的算法，使N！可被素数除以幂

问题描述

是否有一种有效的算法来计算最小整数N，例如N！被p ^ k整除，其中p是一个相对较小的素数，k是一个非常大的整数。换一种说法，

factorial(N) mod p^k == 0

如果给定N和p，我想找出p划分为N！的次数，则可以使用众所周知的公式

k = Sum(floor(N/p^i) for i=1,2,...

我已经用蛮力搜索了较小的k值，但是随着k的增加，这种方法很快就会崩溃，而且似乎没有一种可以推断出较大值的模式。 编辑6/13/2011 通过使用Fiver和Hammar提出的建议，我使用了准二进制搜索来解决问题，但并没有完全按照他们的建议进行。使用上面第二个公式的截断形式，我计算了N的上限，作为k和p的乘积（仅使用第一项）。我使用1作为下限。使用经典的二进制搜索算法，我计算了这两个值之间的中点，并使用第二个公式中的所有中点，将该中点值用作N在第二个公式中计算出k。 如果计算的k太小，我调整下限并重复。太大了，我首先测试看看在中点1处计算出的k是否小于所需的k。如果是这样，则将中点作为最接近的N返回。否则，我调整了高点并重复进行。 如果计算的k相等，则我测试了中点1的值是否等于中点的值。如果是这样，我将高点调整为中点并重复进行。如果中点1小于期望的k，则将中点作为期望的答案返回。 即使k的值很大（10个或更多数字），此方法也可以达到O（n log（n））的速度。     

解决方法

        使用您提到的公式，在给定固定

p

和

N = 1,2...

的情况下，

k

值的序列是不变的。这意味着您可以使用二进制搜索的变体在给定的

k

的情况下找到

N

。 从

N = 1

开始，计算calculate2ѭ。 将

N

加倍，直到

k

大于或等于所需的

k

以获得上限。 对剩余间隔进行二分查找，找到

k

。     ,        好吧，这很有趣。 定义f（i）=（p ^ i-1）/（p-1） 用一种有趣的\“ base \”来写k，其中位置i的值就是这个f（i）。 您可以从最高有效位数到最低有效位数进行此操作。因此，首先找到最大的j，使f（j）<= k。然后计算k / f（j）的商和余数。将商存储为q_j，余数存储为r。现在计算r / f（j-1）的商和余数。将商存储为q_ {j-1}，余数存储为r。现在计算r / f（j-2）的商和余数。等等。 这将生成序列q_j，q_ {j-1}，q_ {j-2}，...，q_1。 （请注意，该序列以1而不是0结尾。）然后计算q_j * p ^ j + q_ {j-1} * p ^（j-1）+ ... q_1 * p。那是你的N。 例如：k = 9，p =3。因此f（i）=（3 ^ i-1）/2。f（1）= 1，f（2）= 4，f（3）= 13。 j与f（j）<= 9等于i = 2与f（2）=4。取9和4的商和余数。这是2的商（这是我们2 \位置的数字）和1的余数。 对于1的余数，找到1 / f（1）的商和余数。商为1，余数为零，所以我们完成了。 所以q_2 = 2，q_1 =1。2 * 3 ^ 2 + 1 * 3 ^ 1 = 21，它是右边的N。 我在纸上有一个解释为什么这样做的原因，但是我不确定如何在文本中进行交流。请注意，f（i）回答了这个问题，“（p ^ i）中有多少个p因子。 ！\”。一旦找到最大的i，j使得j * f（i）小于k，并且意识到自己真正在做的就是找到小于N的最大j * p ^ i，剩下的就从洗。例如，在我们的p = 3例子中，我们得到4 p \'s由1-9的乘积贡献，另外4 p's由10-18的乘积贡献，而另外21由21贡献。前两个只是倍数。的p ^ 2; f（2）= 4告诉我们p ^ 2的每个倍数对乘积贡献另外4个p'。 [更新] 代码总是有助于澄清。将以下perl脚本另存为

foo.pl

，并将其运行为

foo.pl <p> <k>

。请注意，

**

是Perl的幂运算符，ѭ16s计算BigInts（无限精度整数）的商和余数，

use bigint

告诉Perl在各处使用BigInts。

#!/usr/bin/env perl

use warnings;
use strict;
use bigint;

@ARGV == 2
    or die \"Usage: $0 <p> <k>\\n\";

my ($p,$k) = map { Math::BigInt->new($_) } @ARGV;

sub f {
    my $i = shift;
    return ($p ** $i - 1) / ($p - 1);
}

my $j = 0;
while (f($j) <= $k) {
    $j++;
}
$j--;

my $N = 0;

my $r = $k;
while ($r > 0) {
    my $val = f($j);
    my ($q,$new_r) = $r->bdiv($val);
    $N += $q * ($p ** $j);
    $r = $new_r;
    $j--;
}

print \"Result: $N\\n\";

exit 0;

    ,为什么不使用提到的第二个公式尝试二进制搜索答案？ 您只需要考虑N的值，对于该值，p除以N，因为如果不这样，则N！和（N-1）！被p的相同次方除，所以N不能是最小的。     ,        考虑 我=（pn）！ 并忽略除p以外的素数因子。结果看起来像 I = pn * pn-1 * pn-2 * ... * p * 1 I = pn +（n-1）+（n-2）+ ... 2 + 1 I = p（n2 + n）/ 2 因此，我们试图找到最小的n （n2 + n）/ 2> = k 如果我记得二次方程式给我们的话 N = pn，其中n> =（sqrt（1 + 8k）-1）/ 2 （PS。有人知道如何在降价促销中显示基本符号吗？） 编辑： 错了让我看看我是否可以挽救它...