回答編集履歴
1
追記
test
CHANGED
|
@@ -283,3 +283,171 @@
|
|
|
283
283
|
|
|
284
284
|
|
|
285
285
|
軽くプロファイリングした感じだと、一番内側のforループで大半の時間を食っているようです。
|
|
286
|
+
|
|
287
|
+
|
|
288
|
+
|
|
289
|
+
|
|
290
|
+
|
|
291
|
+
### 追記
|
|
292
|
+
|
|
293
|
+
最終的に、numbaのちからで区間篩の方が速くなるところまで持っていきました。numbaのjitコンパイルを通すためにコードをまただいぶ書き換えてあります。
|
|
294
|
+
|
|
295
|
+
```python
|
|
296
|
+
|
|
297
|
+
import timeit
|
|
298
|
+
|
|
299
|
+
import math
|
|
300
|
+
|
|
301
|
+
|
|
302
|
+
|
|
303
|
+
import numba
|
|
304
|
+
|
|
305
|
+
|
|
306
|
+
|
|
307
|
+
#旧：エラトステネスの篩
|
|
308
|
+
|
|
309
|
+
@numba.jit(nopython=True)
|
|
310
|
+
|
|
311
|
+
def eratosthenes_list(n):
|
|
312
|
+
|
|
313
|
+
primes = [True] * n
|
|
314
|
+
|
|
315
|
+
result = [2]
|
|
316
|
+
|
|
317
|
+
for prime in range(3, n, 2):
|
|
318
|
+
|
|
319
|
+
if primes[prime]:
|
|
320
|
+
|
|
321
|
+
result.append(prime)
|
|
322
|
+
|
|
323
|
+
for i in range(prime * prime, n, prime * 2):
|
|
324
|
+
|
|
325
|
+
primes[i] = False
|
|
326
|
+
|
|
327
|
+
return result
|
|
328
|
+
|
|
329
|
+
|
|
330
|
+
|
|
331
|
+
#新：区間篩
|
|
332
|
+
|
|
333
|
+
@numba.jit(nopython=True)
|
|
334
|
+
|
|
335
|
+
def primesList(n):
|
|
336
|
+
|
|
337
|
+
n -= 1
|
|
338
|
+
|
|
339
|
+
SIEVE_MAX = int(math.sqrt(n)) + 1
|
|
340
|
+
|
|
341
|
+
|
|
342
|
+
|
|
343
|
+
sieve = [True] * (SIEVE_MAX)
|
|
344
|
+
|
|
345
|
+
sieve[0] = sieve[1] = False
|
|
346
|
+
|
|
347
|
+
for i in range(2, SIEVE_MAX):
|
|
348
|
+
|
|
349
|
+
for j in range(i * 2, SIEVE_MAX, i):
|
|
350
|
+
|
|
351
|
+
sieve[j] = False
|
|
352
|
+
|
|
353
|
+
|
|
354
|
+
|
|
355
|
+
primes = []
|
|
356
|
+
|
|
357
|
+
for i, x in enumerate(sieve):
|
|
358
|
+
|
|
359
|
+
if x:
|
|
360
|
+
|
|
361
|
+
primes.append(i)
|
|
362
|
+
|
|
363
|
+
|
|
364
|
+
|
|
365
|
+
for i in range(1, (n // (SIEVE_MAX))+1):
|
|
366
|
+
|
|
367
|
+
sieve = [True] * len(sieve)
|
|
368
|
+
|
|
369
|
+
start_index = SIEVE_MAX * i
|
|
370
|
+
|
|
371
|
+
|
|
372
|
+
|
|
373
|
+
for p in primes:
|
|
374
|
+
|
|
375
|
+
if (p * p > start_index + SIEVE_MAX):
|
|
376
|
+
|
|
377
|
+
break
|
|
378
|
+
|
|
379
|
+
|
|
380
|
+
|
|
381
|
+
jj = (p - (start_index % p)) % p
|
|
382
|
+
|
|
383
|
+
for j in range(jj, SIEVE_MAX, p):
|
|
384
|
+
|
|
385
|
+
sieve[j] = False
|
|
386
|
+
|
|
387
|
+
|
|
388
|
+
|
|
389
|
+
tmp = n - start_index + 1
|
|
390
|
+
|
|
391
|
+
if tmp < SIEVE_MAX:
|
|
392
|
+
|
|
393
|
+
x = tmp
|
|
394
|
+
|
|
395
|
+
else:
|
|
396
|
+
|
|
397
|
+
x = SIEVE_MAX
|
|
398
|
+
|
|
399
|
+
for j in range(x):
|
|
400
|
+
|
|
401
|
+
if sieve[j]:
|
|
402
|
+
|
|
403
|
+
primes.append(start_index + j)
|
|
404
|
+
|
|
405
|
+
return primes
|
|
406
|
+
|
|
407
|
+
|
|
408
|
+
|
|
409
|
+
if True:
|
|
410
|
+
|
|
411
|
+
for nn in range(10, 10**3):
|
|
412
|
+
|
|
413
|
+
a = eratosthenes_list(nn)
|
|
414
|
+
|
|
415
|
+
b = primesList(nn)
|
|
416
|
+
|
|
417
|
+
if not a == b:
|
|
418
|
+
|
|
419
|
+
print(nn)
|
|
420
|
+
|
|
421
|
+
print("miss")
|
|
422
|
+
|
|
423
|
+
print(len(a), len(b))
|
|
424
|
+
|
|
425
|
+
print(a)
|
|
426
|
+
|
|
427
|
+
print(b)
|
|
428
|
+
|
|
429
|
+
break
|
|
430
|
+
|
|
431
|
+
else:
|
|
432
|
+
|
|
433
|
+
print("all ok")
|
|
434
|
+
|
|
435
|
+
|
|
436
|
+
|
|
437
|
+
nn = 10**8
|
|
438
|
+
|
|
439
|
+
print(timeit.timeit(lambda: (eratosthenes_list(nn)), number = 1))
|
|
440
|
+
|
|
441
|
+
print(timeit.timeit(lambda: (primesList(nn)), number = 1))
|
|
442
|
+
|
|
443
|
+
""" =>
|
|
444
|
+
|
|
445
|
+
all ok
|
|
446
|
+
|
|
447
|
+
1.2398602159810252
|
|
448
|
+
|
|
449
|
+
1.059988280001562
|
|
450
|
+
|
|
451
|
+
"""
|
|
452
|
+
|
|
453
|
+
```
|