我试图通过去掉子列表长度的复杂公式来进一步优化素数线程中的冠军解决方案。相同子序列的len()太慢,因为len的开销很大,而且生成子序列的开销也很大。这看起来会稍微提高函数的速度,但我还不能去掉除法,尽管我只在条件语句中做了除法。当然,我可以尝试通过对n而不是n*n进行起始标记的优化来简化长度计算。
我将除法/替换为整数除法//以兼容Python 3或
from __future__ import division
如果这个递归公式可以帮助加速numpy解,我也会很感兴趣,但我没有太多使用numpy的经验。
如果您为代码启用了special,那么故事就会完全不同,而且Atkins sieve代码会比这种特殊的切片技术更快。
import cProfile
def rwh_primes1(n):
# http://stackoverflow.com/questions/2068372/fastest-way-to-list-all-primes-below-n-in-python/3035188#3035188
""" Returns a list of primes < n """
sieve = [True] * (n//2)
for i in xrange(3,int(n**0.5)+1,2):
if sieve[i//2]:
sieve[i*i//2::i] = [False] * ((n-i*i-1)//(2*i)+1)
return [2] + [2*i+1 for i in xrange(1,n/2) if sieve[i]]
def primes(n):
# http://stackoverflow.com/questions/2068372/fastest-way-to-list-all-primes-below-n-in-python/3035188#3035188
# recurrence formula for length by amount1 and amount2 Tony Veijalainen 2010
""" Returns a list of primes < n """
sieve = [True] * (n//2)
amount1 = n-10
amount2 = 6
for i in xrange(3,int(n**0.5)+1,2):
if sieve[i//2]:
## can you make recurrence formula for whole reciprocal?
sieve[i*i//2::i] = [False] * (amount1//amount2+1)
amount1-=4*i+4
amount2+=4
return [2] + [2*i+1 for i in xrange(1,n//2) if sieve[i]]
numprimes=1000000
print('Profiling')
cProfile.Profile.bias = 4e-6
for test in (rwh_primes1, primes):
cProfile.run("test(numprimes)")
分析(不同版本之间差别不大)
3 function calls in 0.191 CPU seconds
Ordered by: standard name
ncalls tottime percall cumtime percall filename:lineno(function)
1 0.006 0.006 0.191 0.191 <string>:1(<module>)
1 0.185 0.185 0.185 0.185 myprimes.py:3(rwh_primes1)
1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects}
3 function calls in 0.192 CPU seconds
Ordered by: standard name
ncalls tottime percall cumtime percall filename:lineno(function)
1 0.006 0.006 0.192 0.192 <string>:1(<module>)
1 0.186 0.186 0.186 0.186 myprimes.py:12(primes)
1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects}
有趣的是,通过将限制增加到10**8,并将计时修饰器添加到删除分析的函数中:
rwh_primes1 took 23.670 s
primes took 22.792 s
primesieve took 10.850 s
有趣的是,如果你不生成素数列表,而是返回筛子本身,时间大约是数字列表版本的一半。
发布于 2011-04-13 17:07:19
你可以做一个轮子优化。2和3的倍数不是素数,所以根本不要存储它们。然后,您可以从5开始,通过以2、4、2、4、2、4等步长递增来跳过2和3的倍数。
下面是它的C++代码。希望这能有所帮助。
void sieve23()
{
int lim=sqrt(MAX);
for(int i=5,bit1=0;i<=lim;i+=(bit1?4:2),bit1^=1)
{
if(!isComp[i/3])
{
for(int j=i,bit2=1;;)
{
j+=(bit2?4*i:2*i);
bit2=!bit2;
if(j>=MAX)break;
isComp[j/3]=1;
}
}
}
}
发布于 2011-04-13 17:33:47
如果您可能决定使用C++来提高速度,我将Python sieve移植到了C++。完整的讨论可以在这里找到:Porting optimized Sieve of Eratosthenes from Python to C++。
在使用g++ -O3
和N=100000000编译的英特尔Q6600 Ubuntu 10.10上,这需要415毫秒。
#include <vector>
#include <boost/dynamic_bitset.hpp>
// http://vault.embedded.com/98/9802fe2.htm - integer square root
unsigned short isqrt(unsigned long a) {
unsigned long rem = 0;
unsigned long root = 0;
for (short i = 0; i < 16; i++) {
root <<= 1;
rem = ((rem << 2) + (a >> 30));
a <<= 2;
root++;
if (root <= rem) {
rem -= root;
root++;
} else root--;
}
return static_cast<unsigned short> (root >> 1);
}
// https://stackoverflow.com/questions/2068372/fastest-way-to-list-all-primes-below-n-in-python/3035188#3035188
// https://stackoverflow.com/questions/5293238/porting-optimized-sieve-of-eratosthenes-from-python-to-c/5293492
template <class T>
void primesbelow(T N, std::vector<T> &primes) {
T i, j, k, sievemax, sievemaxroot;
sievemax = N/3;
if ((N % 6) == 2) sievemax++;
sievemaxroot = isqrt(N)/3;
boost::dynamic_bitset<> sieve(sievemax);
sieve.set();
sieve[0] = 0;
for (i = 0; i <= sievemaxroot; i++) {
if (sieve[i]) {
k = (3*i + 1) | 1;
for (j = k*k/3; j < sievemax; j += 2*k) sieve[j] = 0;
for (j = (k*k+4*k-2*k*(i&1))/3; j < sievemax; j += 2*k) sieve[j] = 0;
}
}
primes.push_back(2);
primes.push_back(3);
for (i = 0; i < sievemax; i++) {
if (sieve[i]) primes.push_back((3*i+1)|1);
}
}
https://stackoverflow.com/questions/3285443
复制相似问题