parent
ad67644eeb
commit
63a76b01f8
|
@ -0,0 +1,38 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from operator import mul
|
||||
from functools import reduce
|
||||
|
||||
def factorize(bound):
|
||||
sieve = [[] for i in range(bound)]
|
||||
for n in range(2, bound):
|
||||
if sieve[n] == []:
|
||||
sieve[n] = [n]
|
||||
for i in range(2, n+1):
|
||||
if i*n >= bound:
|
||||
break
|
||||
sieve[i*n] = sieve[i] + sieve[n]
|
||||
return sieve
|
||||
|
||||
def bin(row, col, fnums):
|
||||
num = [f for a in fnums[row-col+1:row+1] for f in a]
|
||||
den = [f for a in fnums[1:col+1] for f in a]
|
||||
for f in set(num):
|
||||
if num.count(f) - den.count(f) > 1:
|
||||
return 0
|
||||
return reduce(mul, num) // reduce(mul, den)
|
||||
|
||||
|
||||
rn = 51
|
||||
fnums = [[1] + i for i in factorize(rn)]
|
||||
t = 0
|
||||
dup = []
|
||||
for row in range(2, rn):
|
||||
for col in range(1, row//2+1):
|
||||
n = bin(row, col, fnums)
|
||||
if n in dup:
|
||||
continue
|
||||
t += n
|
||||
dup.append(n)
|
||||
|
||||
print(t+1)
|
|
@ -0,0 +1,22 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
def primes(n):
|
||||
sieve = [False]*n
|
||||
for p in range(2, n):
|
||||
if sieve[p]:
|
||||
continue
|
||||
for i in range(p**2, n, p):
|
||||
sieve[i] = True
|
||||
return [p for p, t in enumerate(sieve) if not t][2:]
|
||||
|
||||
def rec(n, nums, i, lim):
|
||||
c = 1
|
||||
for j, f in enumerate(nums[i:], start=i):
|
||||
new = f*n
|
||||
if new > lim:
|
||||
break
|
||||
c += rec(new, nums, j, lim)
|
||||
return c
|
||||
|
||||
nums = primes(100)
|
||||
print(rec(1, nums, 0, 10**9))
|
|
@ -0,0 +1,17 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from itertools import product
|
||||
|
||||
peter = [0]*(4*9 + 1)
|
||||
colin = [0]*(6*6 + 1)
|
||||
|
||||
for p in product([1,2,3,4], repeat=9):
|
||||
peter[sum(p)] += 1
|
||||
for p in product([1,2,3,4,5,6], repeat=6):
|
||||
colin[sum(p)] += 1
|
||||
|
||||
c = 0.0
|
||||
for i, p in enumerate(peter):
|
||||
c += p * sum(colin[:i])
|
||||
|
||||
print("{0:.7f}".format(c/(6**6 * 4**9)))
|
|
@ -0,0 +1,23 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
# Brute force
|
||||
from math import sqrt
|
||||
|
||||
def main():
|
||||
low = 10*(int(sqrt(10203040506070809))//10)
|
||||
upp = int(sqrt(19293949596979899))
|
||||
for i in range(low+3, upp+1, 10):
|
||||
if is_num(i**2):
|
||||
return i*10
|
||||
for i in range(low+7, upp+1, 10):
|
||||
if is_num(i**2):
|
||||
return i*10
|
||||
return False
|
||||
|
||||
def is_num(n):
|
||||
s = str(n)
|
||||
for i in range(9):
|
||||
if int(s[2*i]) != i+1: return False
|
||||
return True
|
||||
|
||||
print(main())
|
|
@ -0,0 +1,42 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from math import ceil
|
||||
|
||||
def is_square(n):
|
||||
return int(n**0.5)**2 == n
|
||||
|
||||
def getpw(n, p):
|
||||
pw = 0
|
||||
while not n % p:
|
||||
pw += 1
|
||||
n //= p
|
||||
return pw
|
||||
|
||||
N = 64 * 10**6
|
||||
|
||||
sieve = [1] * N
|
||||
|
||||
for p in range(2, N, 2):
|
||||
sieve[p] = 2
|
||||
|
||||
for p in range(3, int(ceil(N**0.5 - 1)), 2):
|
||||
if sieve[p] > 1:
|
||||
continue
|
||||
sieve[p] = p
|
||||
for i in range(p**2, N, p):
|
||||
sieve[i] = p
|
||||
|
||||
print("Sieved")
|
||||
s = 1
|
||||
|
||||
for n in range(2, N):
|
||||
if sieve[n] == 1:
|
||||
sieve[n] = n
|
||||
p = sieve[n]
|
||||
pw = p ** getpw(n, p)
|
||||
sieve[n] = sieve[n // p] + \
|
||||
pw**2 * sieve[n // pw]
|
||||
if is_square(sieve[n]):
|
||||
s += n
|
||||
|
||||
print(s)
|
|
@ -0,0 +1,40 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
def next(sums, counts, mod=0):
|
||||
ln = len(sums) + 9
|
||||
news = [0] * ln
|
||||
newc = [0] * ln
|
||||
for i, (s, c) in enumerate(zip(sums, counts)):
|
||||
for n in range(10):
|
||||
news[i + n] += 10 * s + c * n
|
||||
newc[i + n] += c
|
||||
return news, newc
|
||||
|
||||
|
||||
left_sums = range(10)
|
||||
left_count = [0] + [1] * 9
|
||||
|
||||
right_sums = range(10)
|
||||
right_count = [1] * 10
|
||||
|
||||
|
||||
length = 47
|
||||
mod = 3**15
|
||||
s = 45
|
||||
|
||||
for nofd in range(2, length+1):
|
||||
partial = sum(l * 10**(nofd//2 + nofd%2) * rc + r * lc
|
||||
for l, r, lc, rc in
|
||||
zip(left_sums, right_sums, left_count, right_count))
|
||||
if not nofd % 2:
|
||||
s += partial
|
||||
s %= mod
|
||||
continue
|
||||
total = sum(lc * rc
|
||||
for lc, rc in zip(left_count, right_count))
|
||||
s += 10 * partial + 45 * 10**(nofd//2) * total
|
||||
s %= mod
|
||||
left_sums, left_count = next(left_sums, left_count)
|
||||
right_sums, right_count = next(right_sums, right_count)
|
||||
|
||||
print(s)
|
|
@ -0,0 +1,42 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from math import ceil
|
||||
|
||||
def reducetwo(n):
|
||||
c = 0
|
||||
while not n%2:
|
||||
n //= 2
|
||||
c += 1
|
||||
return n, c
|
||||
|
||||
def primes(n):
|
||||
sieve = [0] * (n//2)
|
||||
for p in range(3, n, 2):
|
||||
if sieve[p//2] > 0:
|
||||
continue
|
||||
for k in range(1, n//p+1, 2):
|
||||
sieve[(p*k)//2] += p + sieve[k//2] - sieve[(p*k)//2]
|
||||
return sieve
|
||||
|
||||
|
||||
n = 2*10**7
|
||||
k = 15*10**6
|
||||
sieve = primes(n)
|
||||
print("Primed.")
|
||||
s = 0
|
||||
|
||||
for m in range(n, n-k, -1):
|
||||
if not m % 2:
|
||||
m, c = reducetwo(m)
|
||||
s += 2*c
|
||||
s += sieve[m//2]
|
||||
|
||||
print("Calculated.")
|
||||
|
||||
for m in range(2, k+1):
|
||||
if not m % 2:
|
||||
m, c = reducetwo(m)
|
||||
s -= 2*c
|
||||
s -= sieve[m//2]
|
||||
|
||||
print(s)
|
|
@ -0,0 +1,48 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from itertools import combinations
|
||||
from operator import mul
|
||||
from functools import reduce
|
||||
|
||||
def factorize(n, primes):
|
||||
return [p for p in primes if not n%p]
|
||||
|
||||
def totient(factors, n):
|
||||
r = 0
|
||||
for i in range(1, len(factors)+1):
|
||||
for f in combinations(factors, i):
|
||||
r += (-1)**i * (n // reduce(mul, f))
|
||||
return n + r
|
||||
|
||||
def get_primes(n):
|
||||
sieve = [False] * n
|
||||
r = []
|
||||
for p in range(2, n):
|
||||
if sieve[p]:
|
||||
continue
|
||||
r.append(p)
|
||||
for i in range(p**2, n, p):
|
||||
sieve[i] = True
|
||||
return r
|
||||
|
||||
a, b = 15499, 94744
|
||||
primes = get_primes(100)
|
||||
n = 1
|
||||
for i, p in enumerate(primes):
|
||||
n *= p
|
||||
t = totient(primes[:i+1], n)
|
||||
if t*b < (n-1)*a:
|
||||
break
|
||||
else:
|
||||
print("Not found")
|
||||
quit()
|
||||
|
||||
n //= primes[i]
|
||||
factors = primes[:i]
|
||||
|
||||
for i in range(2, 100):
|
||||
t = totient(set(factors + factorize(i, primes)), n*i)
|
||||
if t*b < (n*i - 1) * a:
|
||||
print(n*i)
|
||||
quit()
|
||||
print("Not found")
|
|
@ -0,0 +1,58 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from lib import primegen
|
||||
from itertools import product
|
||||
|
||||
segments = {"th" : (0,1), "mh" : (1,1),
|
||||
"bh" : (2,1), "tl" : (1,0),
|
||||
"bl" : (2,0), "tr" : (1,2),
|
||||
"br" : (2,2)
|
||||
}
|
||||
|
||||
def rootdiff(n):
|
||||
chain = list(rootchain(n))
|
||||
return sum(transdiff(n1, n2)
|
||||
for n1, n2 in zip(chain, chain[1:]))
|
||||
|
||||
def rootchain(n):
|
||||
n = str(n)
|
||||
while len(n) > 1:
|
||||
yield n
|
||||
n = str(sum(map(int, n)))
|
||||
yield n
|
||||
|
||||
def transdiff(n1, n2):
|
||||
return sum(cache[(d1, d2)]
|
||||
for d1, d2 in zip(n1[::-1], n2[::-1]))
|
||||
|
||||
def init_cache(digs):
|
||||
return {(d1, d2) : 2 * len(digs[d1] & digs[d2])
|
||||
for d1, d2 in product(digs, repeat=2)
|
||||
}
|
||||
|
||||
def parsenum(text):
|
||||
return set(s for s, (rw, cl) in segments.items()
|
||||
if str(text[rw] + 3*' ')[cl] != ' '
|
||||
)
|
||||
|
||||
def parse(fname):
|
||||
with open(fname) as f:
|
||||
r = f.read().split('\n')
|
||||
for n in range(10):
|
||||
yield str(n), parsenum(r[3*n : 3*n+3])
|
||||
|
||||
def main(A, B):
|
||||
s = 0
|
||||
gen = primegen(B)
|
||||
for p in gen:
|
||||
if p > A: break
|
||||
s = rootdiff(p)
|
||||
for p in gen:
|
||||
s += rootdiff(p)
|
||||
return s
|
||||
|
||||
digs = {d : segs for d, segs in parse("data/digs.txt")}
|
||||
digs[' '] = set()
|
||||
cache = init_cache(digs)
|
||||
|
||||
print(main(10**7, 2*10**7))
|
|
@ -0,0 +1,42 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from lib import primegen
|
||||
from fractions import Fraction
|
||||
|
||||
def is_prime(n):
|
||||
return n in primes
|
||||
|
||||
def Pprob(n):
|
||||
if is_prime(n):
|
||||
return Fraction(2, 3)
|
||||
return Fraction(1, 3)
|
||||
|
||||
def croak(n):
|
||||
return {'P' : Pprob(n),
|
||||
'N' : 1-Pprob(n)}
|
||||
|
||||
def jump(sqrs, let):
|
||||
r = [0] * len(sqrs)
|
||||
for i in range(1, len(sqrs)-1):
|
||||
prob = croak(i+1)[let] * \
|
||||
Fraction(1, 2) * sqrs[i]
|
||||
r[i+1] += prob
|
||||
r[i-1] += prob
|
||||
r[1] += sqrs[0] * croak(1)[let]
|
||||
r[-2] += sqrs[-1] * croak(len(sqrs))[let]
|
||||
return r
|
||||
|
||||
def start(n):
|
||||
return [Fraction(1, n)] * n
|
||||
|
||||
|
||||
seq = "PPPPNNPPPNPPNPN"
|
||||
N = 500
|
||||
|
||||
primes = set(primegen(N))
|
||||
sqrs = start(N)
|
||||
|
||||
for letter in seq:
|
||||
sqrs = jump(sqrs, letter)
|
||||
|
||||
print(sum(sqrs))
|
|
@ -0,0 +1,16 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
bound = 10**12
|
||||
nums = set()
|
||||
nums.add(1)
|
||||
|
||||
for base in range(2, int(bound**0.5)+1):
|
||||
r = base**2 + base + 1
|
||||
pw = 3
|
||||
while r < bound:
|
||||
if r not in nums:
|
||||
nums.add(r)
|
||||
r += base**pw
|
||||
pw += 1
|
||||
|
||||
print(sum(nums))
|
|
@ -0,0 +1,29 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from lib import primegen
|
||||
from math import log
|
||||
|
||||
def maxpw(n, logN):
|
||||
return int(logN / log(n))
|
||||
|
||||
def M(p, q, logN):
|
||||
return max(q**pw * p**maxpw(p, logN - pw*log(q))
|
||||
for pw in range(1, maxpw(q, logN - log(p))+1)
|
||||
)
|
||||
|
||||
N = 10**7
|
||||
logN = log(N)
|
||||
|
||||
primes = list(primegen(N//2+1))
|
||||
|
||||
s = 0
|
||||
|
||||
for i, p in enumerate(primes):
|
||||
if p > int(N**0.5):
|
||||
break
|
||||
for q in primes[i+1:]:
|
||||
if p*q > N:
|
||||
break
|
||||
s += M(p, q, logN)
|
||||
|
||||
print(s)
|
|
@ -0,0 +1,54 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from itertools import combinations
|
||||
from operator import mul
|
||||
|
||||
def candidates(sieve):
|
||||
n = len(sieve)
|
||||
for i in range(4, n, 4):
|
||||
sieve[i] = True
|
||||
for p in range(3, n, 2):
|
||||
if sieve[p]:
|
||||
continue
|
||||
for i in range(p**2, n, 2*p):
|
||||
sieve[i] = p
|
||||
|
||||
for i in range(2*p**2, n, 2*p**2):
|
||||
sieve[i] = True
|
||||
if not sieve[p-1]:
|
||||
yield p-1
|
||||
|
||||
def is_prime(n):
|
||||
return not sieve[n] and n%2
|
||||
|
||||
def factors(n):
|
||||
r = []
|
||||
while not n%2:
|
||||
r.append(2)
|
||||
n //= 2
|
||||
while sieve[n] != False:
|
||||
r.append(sieve[n])
|
||||
n //= sieve[n]
|
||||
return r + [n]
|
||||
|
||||
def test_alt(n):
|
||||
fctrs = [f for f in factors(n)]
|
||||
for i in range(len(fctrs)//2):
|
||||
for c in combinations(fctrs, i+1):
|
||||
d = reduce(mul, c)
|
||||
if not is_prime(d + n//d):
|
||||
return False
|
||||
return True
|
||||
|
||||
def test(n, fctrs, p=1, i=0, depth=1):
|
||||
if depth > len(fctrs) // 2:
|
||||
return True
|
||||
for j, f in enumerate(fctrs[i:]):
|
||||
if not is_prime(p*f + n//(p*f)) or not test(n, fctrs, p*f, j+1, depth+1):
|
||||
return False
|
||||
return True
|
||||
|
||||
N = 10**8 + 2
|
||||
sieve = [False] * N
|
||||
|
||||
print(sum(c for c in candidates(sieve) if test(c, factors(c))) + 1)
|
|
@ -0,0 +1,10 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from lib import primegen
|
||||
|
||||
def S(p):
|
||||
return (9 * pow(-24, p-2, p)) % p
|
||||
|
||||
N = 10**8
|
||||
|
||||
print(sum(S(p) for p in list(primegen(N))[2:]))
|
|
@ -0,0 +1,44 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from random import randint
|
||||
|
||||
DEPTH = 10
|
||||
|
||||
def is_prime(n):
|
||||
if n == 2:
|
||||
return True
|
||||
if not n % 2 or n < 2:
|
||||
return False
|
||||
for i in range(DEPTH):
|
||||
a = randint(2, n-1)
|
||||
if pow(a, n-1, n) != 1:
|
||||
return False
|
||||
return True
|
||||
|
||||
def harshad_nums(N, n=0, s=0):
|
||||
r = []
|
||||
for i in range(n==0, 10):
|
||||
m = 10*n + i
|
||||
nws = s + i
|
||||
if m > N:
|
||||
break
|
||||
if m % nws:
|
||||
continue
|
||||
r.append((m, nws))
|
||||
r += harshad_nums(N, m, nws)
|
||||
return r
|
||||
|
||||
def main(N):
|
||||
for n, s in harshad_nums(N // 10):
|
||||
if not is_prime(n // s):
|
||||
continue
|
||||
for d in range(1, 10):
|
||||
m = 10*n + d
|
||||
if m >= N:
|
||||
break
|
||||
if is_prime(m):
|
||||
yield m
|
||||
|
||||
N = 10**14
|
||||
|
||||
print(sum(main(N)))
|
|
@ -0,0 +1,21 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from lib import primegen
|
||||
|
||||
def fpw(p, n):
|
||||
pw = 0
|
||||
pwofp = p
|
||||
while pwofp <= n:
|
||||
pw += n // pwofp
|
||||
pwofp *= p
|
||||
return pw
|
||||
|
||||
N = 10**8
|
||||
mod = 10**9 + 9
|
||||
|
||||
k = 1
|
||||
for p in primegen(N+1):
|
||||
k *= pow(p, 2 * fpw(p, N), mod) + 1
|
||||
k %= mod
|
||||
|
||||
print(k)
|
|
@ -0,0 +1,20 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from fractions import Fraction
|
||||
|
||||
def pick(state, depth, numofc, cache={}):
|
||||
if sum(state) == depth:
|
||||
return len(state) - state.count(0)
|
||||
if state in cache:
|
||||
return cache[state]
|
||||
r = 0
|
||||
for c, num in enumerate(state):
|
||||
nws = sorted(state[:c] + (num + 1,) + state[c+1:])
|
||||
nws = tuple(nws)
|
||||
r += Fraction(numofc-num, numofc*len(state) - sum(state)) * \
|
||||
pick(nws, depth, numofc, cache)
|
||||
cache[state] = r
|
||||
return r
|
||||
|
||||
r = pick(7*(0,), 20, 10)
|
||||
print(float(round(r, 9)))
|
|
@ -0,0 +1,48 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from fractions import gcd
|
||||
from itertools import combinations_with_replacement as combs
|
||||
from math import sqrt
|
||||
|
||||
|
||||
def points(a, b, c, d):
|
||||
return (gcds[a][b] + gcds[b][c] + gcds[c][d] + gcds[d][a])//2 + 1
|
||||
|
||||
def is_square(n):
|
||||
return sieve[n]
|
||||
|
||||
def refs(a,b,c,d):
|
||||
r = [ (a,b,c,d),
|
||||
(d,a,b,c),
|
||||
(c,d,a,b),
|
||||
(b,c,d,a),
|
||||
(a,d,c,b),
|
||||
(b,a,d,c),
|
||||
(c,b,a,d),
|
||||
(d,c,b,a)]
|
||||
return list(set(r))
|
||||
|
||||
count = 0
|
||||
m = 100
|
||||
space = range(1, m+1)
|
||||
|
||||
gcds = [None for n in range(m+1)]
|
||||
for a in range(1, m+1):
|
||||
gcds[a] = [0] + [a*b - gcd(a,b) for b in range(1, m+1)]
|
||||
|
||||
sieve = [False] * 2*m**2
|
||||
n = 1
|
||||
while n**2 < len(sieve):
|
||||
sieve[n**2] = True
|
||||
n += 1
|
||||
|
||||
for a,b,c,d in combs(space, 4):
|
||||
perms = [(a,b,c,d), (a,c,b,d), (a,c,d,b)]
|
||||
wins = set()
|
||||
for p in perms:
|
||||
if not is_square(points(*p)):
|
||||
continue
|
||||
wins |= set(refs(*p))
|
||||
count += len(wins)
|
||||
|
||||
print(count)
|
|
@ -0,0 +1,38 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
def getprimes(n):
|
||||
sieve = [True] * n
|
||||
for p in range(3, n, 2):
|
||||
if not sieve[p]:
|
||||
continue
|
||||
for i in range(p**2, n, 2*p):
|
||||
sieve[i] = False
|
||||
return sieve
|
||||
|
||||
def is_prime(n):
|
||||
if n == 2:
|
||||
return True
|
||||
return n % 2 and sieve[n]
|
||||
|
||||
N = 10**8
|
||||
sieve = getprimes(N)
|
||||
primes = [2] + [n for n in range(3, N, 2) if sieve[n]]
|
||||
|
||||
for n in range(2, int(N**0.5)+1):
|
||||
b = n**2
|
||||
if n % 2:
|
||||
b *= 2
|
||||
for i in range(b, N, b):
|
||||
sieve[i] = n
|
||||
|
||||
s = 0
|
||||
|
||||
for p in primes:
|
||||
d = sieve[p+1]
|
||||
k = d + 1
|
||||
while ((p+1) * k**2)/d**2 <= N:
|
||||
if is_prime(((p+1)*k)//d - 1) and is_prime(((p+1)*k**2)//d**2 - 1):
|
||||
s += p + ((p+1)*k)//d - 1 + ((p+1)*k**2)//d**2 - 1
|
||||
k += 1
|
||||
|
||||
print(s)
|
|
@ -0,0 +1,30 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
def modpw(n, p):
|
||||
r = 0
|
||||
while n >= p:
|
||||
r += n // p
|
||||
n //= p
|
||||
return r
|
||||
|
||||
n = 10**8 + 1
|
||||
sieve = [0] * n
|
||||
|
||||
for p in range(2, n):
|
||||
if sieve[p] > 0:
|
||||
continue
|
||||
minN = p
|
||||
frompw = 1
|
||||
topw = 2
|
||||
k = 1
|
||||
while k*p < n:
|
||||
for pw in range(frompw, topw):
|
||||
k *= p
|
||||
for m in range(k, n, k):
|
||||
sieve[m] = max(sieve[m], minN)
|
||||
minN += p
|
||||
frompw = topw
|
||||
topw = modpw(minN, p) + 1
|
||||
|
||||
print(sum(sieve[2:]))
|
||||
|
|
@ -0,0 +1,30 @@
|
|||
_
|
||||
| |
|
||||
|_|
|
||||
|
||||
|
|
||||
|
|
||||
_
|
||||
_|
|
||||
|_
|
||||
_
|
||||
_|
|
||||
_|
|
||||
|
||||
|_|
|
||||
|
|
||||
_
|
||||
|_
|
||||
_|
|
||||
_
|
||||
|_
|
||||
|_|
|
||||
_
|
||||
| |
|
||||
|
|
||||
_
|
||||
|_|
|
||||
|_|
|
||||
_
|
||||
|_|
|
||||
_|
|
|
@ -0,0 +1,68 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
from random import randint
|
||||
from collections import Counter
|
||||
|
||||
FACTOR_SIEVE = []
|
||||
|
||||
def primegen(n):
|
||||
sieve = [False] * n
|
||||
yield 2
|
||||
for p in range(3, n, 2):
|
||||
if sieve[p]:
|
||||
continue
|
||||
for i in range(p**2, n, 2*p):
|
||||
sieve[i] = True
|
||||
yield p
|
||||
|
||||
def iter_factors(n):
|
||||
sieve = init_factors(n)
|
||||
for i in range(n):
|
||||
yield pfactors(i, sieve)
|
||||
|
||||
def init_factors(n):
|
||||
sieve = [False] * n
|
||||
for p in range(3, n, 2):
|
||||
if sieve[p]:
|
||||
continue
|
||||
for i in range(p**2, n, 2*p):
|
||||
sieve[i] = p
|
||||
return sieve
|
||||
|
||||
def pfactors(n, sieve=FACTOR_SIEVE):
|
||||
if n < 1:
|
||||
if not n:
|
||||
return Counter()
|
||||
n = abs(n)
|
||||
fctrs = Counter()
|
||||
while not n % 2:
|
||||
fctrs[2] += 1
|
||||
n //= 2
|
||||
while sieve[n]:
|
||||
fctrs[sieve[n]] += 1
|
||||
n //= sieve[n]
|
||||
if n != 1:
|
||||
fctrs[n] += 1
|
||||
return fctrs
|
||||
|
||||
|
||||
def miller_rabin(n, k=10):
|
||||
if not (n % 2 and n % 3):
|
||||
return False
|
||||
d = n - 1
|
||||
r = 0
|
||||
while not r % 2:
|
||||
d //= 2
|
||||
r += 1
|
||||
for i in range(k):
|
||||
a = randint(2, n-2)
|
||||
x = pow(a, d, n)
|
||||
if x == 1 or x == n - 1:
|
||||
continue
|
||||
for j in range(r-1):
|
||||
x = pow(x, 2, n)
|
||||
if x == n - 1:
|
||||
break
|
||||
else:
|
||||
return False
|
||||
return True
|
Loading…
Reference in New Issue