Solve more problems road to 200.

python/e146.py

@@ -0,0 +1,49 @@
+from lib_prime import is_prime_rabin_miller
+import concurrent.futures
+
+
+def follows_prime_pattern(n):
+    n2 = n * n
+    prev_prime = None
+    for i in [1, 3, 7, 9, 13, 27]:
+        p = n2 + i
+        if not is_prime_rabin_miller(p):
+            return False
+        if prev_prime is not None:
+            for x in range(prev_prime + 2, p, 2):
+                if is_prime_rabin_miller(x):
+                    return False
+        prev_prime = p
+    return True
+
+
+def check_range(n_min, n_max):
+    s = 0
+    for n in range(n_min - (n_min % 10), n_max, 10):
+        if n % 3 == 0:
+            continue
+        if follows_prime_pattern(n):
+            print(n)
+            s += n
+    return s
+
+
+def euler_146():
+    s = 0
+    m = 150_000_000
+    range_size = 500_000
+
+    with concurrent.futures.ProcessPoolExecutor() as executor:
+        futures = []
+        for start in range(1, m, range_size):
+            end = min(start + range_size, m)
+            futures.append(executor.submit(check_range, start, end))
+        for future in concurrent.futures.as_completed(futures):
+            s += future.result()
+    return s
+
+
+if __name__ == "__main__":
+    solution = euler_146()
+    print("e146.py: " + str(solution))
+    assert solution == 676333270

python/e148.py

@@ -0,0 +1,20 @@
+def count(i, s, c, depth):
+    if depth == 0:
+        return s + i, c + 1
+    for i in range(i, i * 8, i):
+        s, c = count(i, s, c, depth - 1)
+        if c == 10**9:
+            break
+    return s, c
+
+
+def euler_148():
+    s, c = count(1, 0, 0, 11)
+    assert c == 10**9
+    return s
+
+
+if __name__ == "__main__":
+    solution = euler_148()
+    print("e148.py: " + str(solution))
+    assert solution == 2129970655314432

python/e164.py

@@ -0,0 +1,25 @@
+from functools import lru_cache
+
+
+@lru_cache
+def count(last, left):
+    if left == 0:
+        return 1
+    c = 0
+    for d in range(10):
+        if sum(last) + d <= 9:
+            c += count((last[-1], d), left - 1)
+    return c
+
+
+def euler_164():
+    c = 0
+    for d in range(1, 10):
+        c += count((d, ), 20-1)
+    return c
+
+
+if __name__ == "__main__":
+    solution = euler_164()
+    print("e164.py: " + str(solution))
+    assert solution == 378158756814587

python/e188.py

@@ -0,0 +1,38 @@
+from lib_prime import prime_factors
+
+
+def phi(n):
+    # Calculate phi using Euler's totient function
+    c = n
+    fs = set(prime_factors(n))
+    for f in fs:
+        # same as c *= (1 - 1/f)
+        c *= (f - 1)
+        c //= f
+    return c
+
+
+def tetmod(a, k, m):
+    # https://stackoverflow.com/questions/30713648/how-to-compute-ab-mod-m
+    assert k > 0
+    if k == 1:
+        return a
+    if m == 1:
+        if a % 2 == 0:
+            return 0
+        else:
+            return 1
+    phi_m = phi(m)
+    return pow(a, phi_m + tetmod(a, k - 1, phi_m), m)
+
+
+def euler_188():
+    assert phi(100) == 40
+    assert tetmod(14, 2016, 100) == 36
+    return tetmod(1777, 1855, 10**8)
+
+
+if __name__ == "__main__":
+    solution = euler_188()
+    print("e188.py: " + str(solution))
+    assert solution == 95962097

python/e195.py

@@ -0,0 +1,50 @@
+import math
+
+
+def inradius(a, b, c):
+    s = (a + b + c) / 2
+    area = (s * (s - a) * (s - b) * (s - c))
+    inradius = area / s**2
+    return inradius
+
+
+def euler_195():
+    count = 0
+    limit = 1053779
+    limit2 = limit * limit
+    for m in range(0, limit2):
+        for n in range(1, m // 2 + 1):
+            if math.gcd(m, n) == 1:
+                # Generator function from wiki for 60 degree eisenstein
+                # triangles.
+                a = m**2 - m*n + n**2
+                b = 2*m*n - n**2
+                c = m**2 - n**2
+                gcd = math.gcd(a, b, c)
+                a, b, c = a // gcd, b // gcd, c // gcd
+
+                # if gcd == 3 and inradius(a, b, c) > limit2:
+                #     break
+                if inradius(a / 3, b / 3, c / 3) > limit2:
+                    if n == 1:
+                        return count
+                    break
+
+                if a == b and b == c:
+                    continue
+
+                k = 1
+                while True:
+                    ir = inradius(a*k, b*k, c*k)
+                    if ir > limit2:
+                        break
+                    count += 1
+                    k += 1
+
+    return count
+
+
+if __name__ == "__main__":
+    solution = euler_195()
+    print("e195.py: " + str(solution))
+    assert solution == 75085391 

python/lib_prime.py

@@ -102,7 +102,6 @@ def prime_nth(n):
     is 1 indexed, i.e. n = 1 -> 2, n = 2 -> 3, etc.
 
     :param n:
-
     """
     if n == 1:
         return 2