在Python中,我有两种查找素数的算法。每个循环的内部循环似乎执行了相同的次数,并且同样简单。但是,其中一个的花费是另一个的十倍。我的问题是:
为什么?这是可以被优化掉的Python怪癖吗(如何?),或者我还缺少其他东西吗?
我要解决的问题基本上来自http://www.spoj.pl/problems/PRIME1/。在我的情况下,我有N = 10 9,delta = 10 5,我想找到N- delta和delta之间的所有素数。我还拥有smallprimes一个小于或等于N的平方根的所有素数的预制列表。第一个算法非常简单:
smallprimes
def range_f1(lo, hi, smallprimes): """Finds all primes p with lo <= p <= hi. smallprimes is the sorted list of all primes up to (at least) square root of hi. hi & lo might be large, but hi-lo+1 miust fit into a long.""" primes =[] for i in xrange(hi-lo+1): n = lo + i isprime = True for p in smallprimes: if n % p == 0: isprime = False break if isprime: primes.append(n) return primes
通话range_f1(N-delta,N,smallprimes)时间很长(约10秒)。内循环称为195170次。我也有此算法的一个版本,该版本用列表理解替换了循环(这是我实际用于性能分析的列表;请参阅问题的结尾),但这并不快。
range_f1(N-delta,N,smallprimes)
第二个版本是Eratosthenes筛子的一个丑陋的实现:
def range_sieve(lo, hi, smallprimes): """Parameters as before""" # So ugly! But SO FAST!! How?? delta = hi-lo+1 iamprime = [True] * delta # iamprime[i] says whether lo + i is prime if lo<= 1: iamprime[1 - lo] = False def sillyfun(): # For profiling & speed comparison pass for p in smallprimes: rem = lo % p if rem == 0: iamprime[0] = False for i in xrange(p - rem, delta, p): iamprime[i] = False sillyfun() if p >= lo and p <= hi: iamprime[p - lo] = True return [p + lo for (p, ami) in enumerate(iamprime) if ami]
这大约快10倍,耗时不到2秒。但是,内部循环(sillyfun())被调用259982次,比上一种情况要多。我无所适从地解释了为什么这样做很快。
我认为可能是因为第一种算法的内部循环包含模块化算法,而第二种算法只有一个赋值。但是,以下内容似乎暗示赋值并不比模数运算快:
>>> from timeit import timeit >>> timeit("10**9 % 2341234") 0.23445186469234613 >>> timeit("a[5000]=False", "a = [True] * 10000") 0.47924750212666822
这是我实际使用的第一个算法的(可读性较低)版本:
def range_f2(lo, hi, smallprimes): primes =[] for i in xrange(hi-lo+1): n = lo + i try: (1 for p in smallprimes if n % p ==0).next() except StopIteration: primes.append(n) return primes
这是为range_f2()调用事件探查器的结果(请注意,计算生成表达式的时间数):
>>> from cProfile import run as prof >>> prof("range_f2(N-delta,N,sp)") 200005 function calls in 13.866 CPU seconds Ordered by: standard name ncalls tottime percall cumtime percall filename:lineno(function) 1 0.000 0.000 13.866 13.866 <string>:1(<module>) 195170 12.632 0.000 12.632 0.000 prime1.py:101(<genexpr>) 1 1.224 1.224 13.865 13.865 prime1.py:90(range_f2) 4832 0.009 0.000 0.009 0.000 {method 'append' of 'list' objects} 1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects}
这是range_sieve()的分析器结果:
>>> prof("range_sieve(N-delta,N,sp)") 259985 function calls in 1.370 CPU seconds Ordered by: standard name ncalls tottime percall cumtime percall filename:lineno(function) 1 0.003 0.003 1.370 1.370 <string>:1(<module>) 1 0.877 0.877 1.367 1.367 prime1.py:39(range_sieve) 259982 0.490 0.000 0.490 0.000 prime1.py:51(sillyfun) 1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects}
UPDATE 根据流行的需求,第一部分代码的概要分析数据。没有关于内部循环执行多少次的数据,但很显然它与第三个相同。
>>> prof("range_f1(N-delta,N,sp)") 4835 function calls in 14.184 CPU seconds Ordered by: standard name ncalls tottime percall cumtime percall filename:lineno(function) 1 0.000 0.000 14.184 14.184 <string>:1(<module>) 1 14.162 14.162 14.183 14.183 prime1.py:69(range_f1) 4832 0.021 0.000 0.021 0.000 {method 'append' of 'list' objects} 1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects}
区别在于算法上的区别。在第一个版本中,您会使用所有小素数来测试每个候选者-当小素数超过上限时,您不会停下来candidate ** 0.5并不重要,因为小素数range(10**9 - 10**5 ,10**9)是否具有良好的上限,但如果范围相对于上限要大得多。对于复合材料,这不会产生太大的成本,因为它们中的大多数都有至少一个小除数。但是对于素数,您必须一直到N**0.5。10**5/log(10**9)在该时间间隔中大约有素数,每个素数均按素数进行试验划分10**4.5/log(10**4.5),因此可以1.47*10**7对素数进行试验划分。
candidate ** 0.5
range(10**9 - 10**5 ,10**9)
N**0.5
10**5/log(10**9)
10**4.5/log(10**4.5)
1.47*10**7
另一方面,使用筛子时,您只能将合成物舍掉,每个合成物被除以质数的次数会被除掉多次,而质数却不会被除掉(因此质数是免费的)。的素数除数n受(的整数倍)log(n)(这是一个粗略的上限,通常会大大高估)的边界,因此得出10**5*log(10**9)(约一个小常数)相交的上限2*10**6。除此之外,相除可能要比除法少(不了解Python,它适用于C数组)。这样您就可以减少筛子的工作量。
n
log(n)
10**5*log(10**9)
2*10**6
编辑:实际收集的号码10**9-10**5来10**9。
10**9-10**5
10**9
Ticks: 259987 4832 Divisions: 20353799 4832
筛分仅进行259987次穿越,您会看到上方的粗略上限被高估了很多。审判部门需要超过2000万个部门,其中素数为16433632(x/log x低估了素数x = 10**4.5约10%),该范围内的3297个数字使用了3435183,最小素因数大于n**(1/3)。
x/log x
x = 10**4.5
n**(1/3)