euler.py

#! /usr/bin/python3
# vim: set tabstop=4 shiftwidth=4 expandtab :
# Generated by letscode

"""Project Euler solutions"""
import string
import math
import itertools
import collections
import heapq
from prime import prime_factors_list, sieve, nb_divisors, yield_primes
from prime import primes_up_to, nb_prime_divisors, totient, divisors_sieve
from prime import is_prime, prime_divisors_sieve, mult
from functions import ceil, gcd, modinv
from functions import fibo, lcmm, yield_pythagorean_triples_of_peri
from functions import Tn, Pn, Hn, isPn, champernowne_digit
from poker import Hand
import os
from timeit import default_timer as timer


resource_folder = ''


def euler1(lim=1000):
    """Solution for problem 1."""
    # could use sum formula here
    return sum(i for i in range(lim) if i % 3 == 0 or i % 5 == 0)


def euler2(lim=4000000):
    """Solution for problem 2."""
    s = 0
    for f in fibo(1, 2):
        if f > lim:
            return s
        if f % 2 == 0:
            s += f


def euler3(n=600851475143):
    """Solution for problem 3."""
    return prime_factors_list(n)[-1]


def euler4(l=3):
    """Solution for problem 4."""
    # simple optimisation would be an early break
    return max(n for n in (i * j for i, j in itertools.combinations(range(10 ** (l - 1), 10 ** l), 2)) if str(n) == str(n)[::-1])


def euler5(lim=20):
    """Solution for problem 5."""
    return lcmm(*range(1, lim + 1))


def euler6(lim=100):
    """Solution for problem 6."""
    # could use sum formula here
    numbers = range(1, lim + 1)
    sum_ = sum(numbers)
    return sum_ * sum_ - sum(i * i for i in numbers)


def nth(iterable, n, default=None):
    """Returns the nth item or a default value.
    From http://stackoverflow.com/questions/12007820/better-ways-to-get-nth-element-from-an-unsubscriptable-iterable ."""
    return next(itertools.islice(iterable, n, None), default)


def euler7(n=10001):
    """Solution for problem 7."""
    return nth(yield_primes(), n - 1)


def euler8(l=13):
    """Solution for problem 8."""
    # Too lazy for optimisations
    n = '7316717653133062491922511967442657474235534919493496983520312774506326239578318016984801869478851843858615607891129494954595017379583319528532088055111254069874715852386305071569329096329522744304355766896648950445244523161731856403098711121722383113622298934233803081353362766142828064444866452387493035890729629049156044077239071381051585930796086670172427121883998797908792274921901699720888093776657273330010533678812202354218097512545405947522435258490771167055601360483958644670632441572215539753697817977846174064955149290862569321978468622482839722413756570560574902614079729686524145351004748216637048440319989000889524345065854122758866688116427171479924442928230863465674813919123162824586178664583591245665294765456828489128831426076900422421902267105562632111110937054421750694165896040807198403850962455444362981230987879927244284909188845801561660979191338754992005240636899125607176060588611646710940507754100225698315520005593572972571636269561882670428252483600823257530420752963450'
    return max(mult(int(c) for c in s) for s in (n[i:i + l] for i in range(len(n))))


def euler9(p=1000):
    """Solution for problem 9."""
    for a, b, c in yield_pythagorean_triples_of_peri(p):
        return a * b * c


def euler10(lim=2000000):
    """Solution for problem 10."""
    return sum(primes_up_to(lim))


def euler11_():
    """Solution for problem 11."""
    pass


def euler12(nb_div=500):
    """Solution for problem 12."""
    # T(n) = n(n+1)/2
    # T(2k + 0) = 2k*(2k+1) / 2 = k*(2k+1)
    #    and gcd(k, 2*k+1) = gcd(k, 1) = 1
    #    Thus, nb_div(T(2k)) = nb_div(k) * nb_div(2k+1)
    # T(2k + 1) = (2k+1)*(2k+2) / 2 = (2k+1)(k+1)
    #    and gcd(2k+1, k+1) = gcd(k, k+1) = gcd(k, 1) = 1
    #    Thus, nb_div(T(2k + 1)) = nb_div(k+1) * nb_div(2k+1)
    nb_div_k, nb_div_kp1 = None, 1
    for k in itertools.count(1):
        nb_div_k, nb_div_kp1 = nb_div_kp1, nb_divisors(k + 1)
        c = 2 * k + 1
        nb_div_c = nb_divisors(c)
        if nb_div_k * nb_div_c >= nb_div:
            return c * k
        if nb_div_kp1 * nb_div_c >= nb_div:
            return c * (k + 1)


def euler13_():
    """Solution for problem 13."""
    pass


def euler14(lim=1000000):
    """Solution for problem 14."""
    collatz = {1: 1}
    for i in range(2, lim):
        if i not in collatz:
            chain = []
            while i not in collatz:
                chain.append(i)
                i = (3 * i + 1) if i % 2 else (i // 2)
            stop = collatz[i] + 1
            for idx, val in enumerate(reversed(chain)):
                collatz[val] = stop + idx
        assert i in collatz
    return max((v, k) for k, v in collatz.items() if k <= lim)[1]


def euler15(col=20, row=20):
    """Solution for problem 15."""
    nb_routes = [[0] * (row + 1) for i in range(col + 1)]
    nb_routes[0][0] = 1
    for i in range(col + 1):
        for j in range(row + 1):
            if i:
                nb_routes[i][j] += nb_routes[i - 1][j]
            if j:
                nb_routes[i][j] += nb_routes[i][j - 1]
    return nb_routes[-1][-1]


def sum_digit(n):
    return sum(int(c) for c in str(n))


def euler16(n=1000):
    """Solution for problem 16."""
    return sum_digit(2 ** n)


def euler17(lim=1000):
    """Solution for problem 17."""
    # Optimisation could be done by grouping similar number
    # prefixes/suffixes and use multiplications instead of
    # repeated sums but this is fast enough
    small_numbers_str = {
        1: "one", 2: "two", 3: "three", 4: "four", 5: "five", 6: "six",
        7: "seven", 8: "eight", 9: "nine", 10: "ten", 11: "eleven",
        12: "twelve", 13: "thirteen", 14: "fourteen", 15: "fifteen",
        16: "sixteen", 17: "seventeen", 18: "eighteen", 19: "nineteen"
    }
    ten_multiples_str = {
        2: "twenty", 3: "thirty", 4: "forty", 5: "fifty",
        6: "sixty", 7: "seventy", 8: "eighty", 9: "ninety"
    }
    small_numbers_len = {k: len(v) for k, v in small_numbers_str.items()}
    ten_multiples_len = {k: len(v) for k, v in ten_multiples_str.items()}
    and_len = len("and")
    hundred_len = len("hundred")
    thousand_len = len("thousand")

    def get_len_for_number(n):
        n, units = divmod(n, 10)
        n, tens = divmod(n, 10)
        thousands, hundreds = divmod(n, 10)
        total = 0
        if thousands:
            total += small_numbers_len[thousands] + thousand_len
        if hundreds:
            total += small_numbers_len[hundreds] + hundred_len
            if tens or units:
                total += and_len
        if tens > 1:
            total += ten_multiples_len[tens]
        else:
            units += 10 * tens  # move tens in units
        if units:
            total += small_numbers_len[units]
        return total
    return sum(get_len_for_number(i) for i in range(1, lim + 1))


def euler18_():
    """Solution for problem 18."""
    pass


def euler19():
    """Solution for problem 19."""
    day = 2  # Tuesday 1 Jan 1901
    nb_days = 7
    freq = [0] * nb_days
    months = [nb % nb_days for nb in (31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31)]
    for year in range(1901, 2001):
        for i, month in enumerate(months):
            freq[day] += 1
            day = (day + month + (1 if (i == 1 and (year % 4 == 0 and (year % 100 != 0 or year % 400 == 0))) else 0)) % nb_days
    return freq[0]


def euler20(n=100):
    """Solution for problem 20."""
    return sum_digit(math.factorial(n))


def euler21(lim=10000):
    """Solution for problem 21."""
    sum_div = [sum(l) for l in divisors_sieve(lim)]
    return sum(
        s + i
        for i, s in enumerate(sum_div)
        if s < i and sum_div[s] == i)


def euler22(f='p022_names.txt'):
    """Solution for problem 22."""
    with open(os.path.join(resource_folder, f)) as file_:
        return sum((i + 1) * sum(1 + ord(c) - ord('A') for c in name)
                   for i, name in enumerate(sorted(''.join(file_.readlines()).replace('"', '').split(','))))


def euler23(lim=28123):
    """Solution for problem 23."""
    abun = [i > 0 and sum(l) > i for i, l in enumerate(divisors_sieve(lim))]
    return sum(
        n for n in range(lim)
        if not any(i for i in range(1 + n // 2) if abun[i] and abun[n - i]))


def euler24():
    """Solution for problem 24."""
    return int(''.join(nth(itertools.permutations(string.digits), 1000000 - 1)))


def euler25(nb_digits=1000):
    """Solution for problem 25."""
    lim = 10 ** (nb_digits - 1)
    for i, f in enumerate(fibo()):
        if f > lim:
            return 1 + i


def length_recur_cycle(a, b, base=10):
    remainders = {}
    for i in itertools.count():
        if a == 0:
            return 0
        if a in remainders:
            return i - remainders[a]
        remainders[a] = i
        a = (a * base) % b


def euler26(lim=1000):
    """Solution for problem 26."""
    return max(range(2, lim), key=lambda i: length_recur_cycle(1, i))


def euler27(lim=1000):
    """Solution for problem 27."""
    # P(0) = b must be prime
    # P(1) = 1 + a + b must be prime so 1+a+b >= 2 => a => 1-b
    maxa, maxb, maxn = None, None, 0
    for b in primes_up_to(lim):
        for a in range(1 - b, lim + 1):
            assert -lim <= a <= lim
            n = 0
            while is_prime(n * n + a * n + b):
                n += 1
            if n > maxn:
                maxa, maxb, maxn = a, b, n
    return maxa * maxb


def euler28(n=1001):
    """Solution for problem 28."""
    # For a level of size s (s>1), corners are :
    # s*s, s*s-s+1, s*s-2s+2, s*s-3s+3
    # Their sum is 4*s*s - 6*s + 6  (3 <= s <= n)
    # Writing s = 2*j + 3, sum of corners is
    # 16*j*j + 36j + 24
    # This sum could be computed in constant time but this is
    # fast enough.
    return 1 + sum(16 * j * j + 36 * j + 24 for j in range(n // 2))


def euler29(lima=100, limb=100):
    """Solution for problem 29."""
    n = set()
    for a in range(2, lima + 1):
        p = a
        for _ in range(2, limb + 1):
            p *= a
            n.add(p)
    return len(n)


def euler30():
    """Solution for problem 30."""
    # sum(fifth(digits)) = n
    # nb_digits * fifth(0) < sum(fifth(digits)) <= nb_digits * fifth(base - 1)
    # base ^ (nb_digits - 1) <= n < base ^ nb_digits
    # => base ^ (nb_digits - 1) <= nb_digits * fifth(base - 1)
    # with base = 10, nb_digits <= 6
    fifth = {c: int(c) ** 5 for c in string.digits}
    dict_sum = {}
    for nb_dig in range(2, 6 + 1):
        for l in itertools.combinations_with_replacement('123456789', nb_dig):
            dict_sum.setdefault(sum(fifth[c] for c in l), []).append(list(l))
    return sum(n for n, l in dict_sum.items() if [c for c in sorted(str(n)) if c != '0'] in l)


def euler31(obj=200, coins=None):
    """Solution for problem 31."""
    if coins is None:
        coins = [1, 2, 5, 10, 20, 50, 100, 200]
    nb_ways = [1] + [0] * obj
    for c in coins:
        for v in range(obj + 1 - c):
            nb_ways[v + c] += nb_ways[v]
    return nb_ways[-1]


def euler32():
    """Solution for problem 32."""
    # nb_dig(a*b) = nb_dig(a) + nb_dig(b) or nb_dig(a) + nb_dig(b) - 1
    # we want : a * b = c with 0 < a < b < c
    # and we must have : nb_dig(a) + nb_dig(b) + nb_dig(c) = 9
    # i + j + k = 9 and (i + j = k or i + j - 1 = k)
    # => 2i + 2j = 9 (impossible) or 2i + 2j = 10
    # => i + j = 10 => (1, 4, 4), (2, 3, 4).
    return sum({a * b
                for a in range(2, 98)
                for b in range(1234 if a <= 10 else 123, 10000 // a + 1)
                if ''.join(sorted(list(str(a) + str(b) + str(a * b)))) == '123456789'
                })


def euler33():
    """Solution for problem 33."""
    # Different options reducing to b/c (with 0 < b < c < 10 and 0 <= a < 10) are :
    # ab / ac (possible only if a=0 or b=c : not interesting)
    # ba / ca (possible only if a=0 or b=c : not interesting)
    # ab / ca (possible iif 9bc = a (10c - b))
    # ba / ac (possible iif 9bc = a (10b - c))
    # Also, top / bottom = b/c <=> top * c = bottom * b
    t, b = [mult(lst)
            for lst in zip(*[(b, c)
                             for b, c in itertools.combinations(range(1, 10), 2)
                             for a in range(1, 10)
                             for top, bot in [(10 * a + b, 10 * c + a), (10 * b + a, 10 * a + c)]
                             if top * c == bot * b])]
    return b // gcd(t, b)


def euler34():
    """Solution for problem 34."""
    # sum(fact(digits)) = n
    # nb_digits * fact(0) < sum(fact(digits)) <= nb_digits * fact(base - 1)
    # base ^ (nb_digits - 1) <= n < base ^ nb_digits
    # => base ^ (nb_digits - 1) <= nb_digits * fact(base - 1)
    # with base = 10, we have nb_digits <= 7
    fact = {c: math.factorial(int(c)) for c in string.digits}
    dict_sum = {}
    for nb_dig in range(2, 7 + 1):
        for l in itertools.combinations_with_replacement(string.digits, nb_dig):
            dict_sum.setdefault(sum(fact[c] for c in l), []).append(list(l))
    return sum(n for n, l in dict_sum.items() if sorted(str(n)) in l)


def euler35(nb_dig_max=6):
    # permutations of 2 digits or more must contain only 1, 3, 7, 9
    count = 4  # counting 2, 3, 5 and 7
    final_numbers = {'1', '3', '7', '9'}
    for l in range(2, nb_dig_max + 1):
        for p in itertools.product(final_numbers, repeat=l):
            p_int = int(''.join(p))
            perm = {int(''.join(p[i:] + p[:i])) for i in range(len(p))}
            if p_int == min(perm) and all(is_prime(n) for n in perm):
                count += len(perm)
    return count


def euler36():
    """Solution for problem 36."""
    # 999999
    sol = []
    for i in range(1000):
        s = str(i)
        for beg in (-1, -2):
            n = int(s + s[beg::-1])
            s2 = bin(n)[2:]
            if s2 == s2[::-1]:
                sol.append(n)
    return sum(sol)


def euler37():
    """Solution for problem 37."""
    left, sol = {0}, []
    while left:
        left = {n for n in (10 * l + i for i in range(10) for l in left) if is_prime(n)}
        sol.extend(n for n in left if n > 10 and all(is_prime(n % (10 ** pow)) for pow in range(1, len(str(n)))))
    return sum(sol)


def euler38():
    """Solution for problem 38."""
    # '123456789' will be decomposed in at least two elements,
    # the smallest being 4 at most characters long
    sol = 0
    digits = {str(d) for d in range(1, 10)}
    for n in range(10000):
        s = ""
        for i in itertools.count(1):
            s += str(n * i)
            if len(s) >= len(digits):
                if len(s) == len(digits) and set(s) == digits:
                    sol = max(sol, int(s))
                break
    return sol


def euler39(lim=1000):
    """Solution for problem 39."""
    return max(range(1, lim + 1), key=lambda p: len(list(yield_pythagorean_triples_of_peri(p))))


def euler40():
    """Solution for problem 40."""
    return mult(champernowne_digit(pow(10, i) - 1) for i in range(6))


def euler41():
    """Solution for problem 41."""
    # sum(i, i=1..n) is is not divisible by 3 only if n = 1, 4, 7 or bigger than 9
    for nb_dig in (7, 4, 1):
        for l in itertools.permutations(str(d) for d in range(nb_dig, 0, -1)):
            n = int(''.join(l))
            if is_prime(n):
                return n


def euler42_():
    """Solution for problem 42."""
    pass


def euler43():
    """Solution for problem 43."""
    # Could be optimised by using a constructive solution starting from the end
    div = [(17, 7), (13, 6), (11, 5), (7, 4), (5, 3), (3, 2), (2, 1)]
    return sum(int(''.join(p))
               for p in itertools.permutations(string.digits)
               if all(int(''.join(p[i:i + 3])) % d == 0 for d, i in div))


def euler44():
    """Solution for problem 44."""
    # We look for minimal value of P(m) = P(k) - P(j)
    # Such that P(k) + P(j) is pentagonal.
    # We must have j < m < k.
    for m in itertools.count():
        pm = Pn(m)
        for j in range(1, m):
            pj = Pn(j)
            pk = pm + pj
            if isPn(pk) and isPn(pj + pk):
                return pm


def euler45(nb_fact=4):
    """Solution for problem 45."""
    t_gen = (Tn(n) for n in itertools.count())
    p_gen = (Pn(n) for n in itertools.count())
    h_gen = (Hn(n) for n in itertools.count())
    t, p, h = next(t_gen), next(p_gen), next(h_gen)
    while True:
        if p > t:
            t = next(t_gen)
        elif t > p or h > p:
            p = next(p_gen)
        elif p > h:
            h = next(h_gen)
        elif t in [0, 1, 40755]:
            t, p, h = next(t_gen), next(p_gen), next(h_gen)
        else:
            return t


def euler46():
    """Solution for problem 46."""
    for i in itertools.count(9, 2):
        if not is_prime(i) and not any(is_prime(i - 2 * n * n) for n in range(1, ceil(math.sqrt(i // 2)))):
            return i


def euler47(nb_fact=4):
    """Solution for problem 47."""
    cand = []
    for i in itertools.count(2):
        if nb_prime_divisors(i) == nb_fact:
            cand.append(i)
            if len(cand) == nb_fact:
                return cand[0]
        else:
            cand = []


def euler48(n=1000, nb_dig=10):
    """Solution for problem 48."""
    mod = 10 ** nb_dig
    return sum(pow(i, i, mod) for i in range(1, n + 1)) % mod


def sorted_number(n):
    # Reversed to keep 0 instead of` discarding them
    return int(''.join(sorted(str(n), reverse=True)))


def euler49(nb_digit=4):
    """Solution for problem 49."""
    low = 10 ** (nb_digit - 1)
    high = 10 ** nb_digit - 1
    prime = sieve(high)
    prime_perm = {}
    for i in range(low, high + 1):
        if prime[i]:
            prime_perm.setdefault(sorted_number(i), []).append(i)
    for perms in prime_perm.values():
        # could iterate only on a,b but not a bottleneck
        for a, b, c in itertools.combinations(perms, 3):
            assert c > b > a
            if b - a == c - b and a != 1487:
                return int(str(a) + str(b) + str(c))


def euler50(lim=1000000):
    """Solution for problem 50."""
    primes = sieve(lim)
    list_primes = [i for i, p in enumerate(primes) if p]
    max_len, max_sum = 0, 0
    for i in range(len(list_primes)):
        for j in range(i + max_len + 1, len(list_primes)):
            s = sum(list_primes[i:j])  # could use sum array here
            if s > lim:
                break
            elif primes[s]:
                assert j - i > max_len
                max_len, max_sum = j - i, s
    return max_sum


def euler51(nb=8):
    """Solution for problem 51."""
    # Iterating on batches of number of same length to be able to use a sieve.
    # For each prime, looking for digits one can change, then the possible
    # combinations of these digits. For each combination, we compute the value
    # of the corresponding mask and check primality.
    base = 10  # Not so much of a variable due to various optimisations
    lim = 1 + base - nb  # replaced digit must be smaller than lim

    # optimisation : precomputing powers of 10
    pows = []
    # optimisation : the mask (and so its length) must be divisible by 3 for
    # big families, otherwise, we'd have too many multiple of 3.
    len_mask_step = 3 if nb >= 8 else 1
    # optimisation : last number can't be changed for families of 4 or more
    last_nb_can_chg = nb <= 4

    low = 1
    while True:
        high = base * low
        primes = sieve(high)
        pows.append(low)
        for n in range(low, high):
            if primes[n]:
                n_str = [int(c) for c in reversed(str(n))]
                for val_replaced in range(lim):
                    pos = [i for i, c in enumerate(n_str)
                           if c == val_replaced and (i or last_nb_can_chg)]
                    for len_mask in range(len_mask_step, 1 + len(pos), len_mask_step):
                        for pos_mask in itertools.combinations(pos, len_mask):
                            val_mask = sum(pows[i] for i in pos_mask)
                            if 1 + sum(1 for i in range(1, base - val_replaced) if primes[n + i * val_mask]) >= nb:
                                return n
        low = high


def euler52(lim=6):
    """Solution for problem 52."""
    for x in itertools.count(1):
        digits = sorted_number(x)
        if all(digits == sorted_number(i * x) for i in range(2, lim + 1)):
            return x


def euler53(n_max=100, lim=1000000):
    """Solution for problem 53."""
    # C(n, r) = n! / (r! (n-r)!) for r <= n
    #        = prod(n + 1 - i / i, i =1..r)
    def C(n, r):
        r = min(r, n - r)  # symmetry
        return mult(n - i for i in range(r)) // mult(i + 1 for i in range(r))
    s = 0
    for n in range(n_max + 1):
        for r in range(n // 2 + 1):  # stop halfway - use symmetry
            if C(n, r) > lim:
                s += 1 if 2 * r == n else 2
    return s


def euler54(f='p054_poker.txt'):
    """Solution for problem 54."""
    ret = 0
    with open(os.path.join(resource_folder, f)) as file_:
        for l in file_:
            hand1, hand2 = Hand.from_string(l)
            ret += hand1 > hand2
    return ret


def euler55_():
    """Solution for problem 55."""
    pass


def euler56(lim=100):
    """Solution for problem 56."""
    int_ = {str(d): d for d in range(10)}
    maxi = 1
    for a in range(2, lim + 1):
        p = a
        for b in range(2, lim + 1):
            p *= a
            maxi = max(maxi, sum(int_[d] for d in str(p)))
    return maxi


def euler57(nb_exp=1000):
    """Solution for problem 57."""
    # If expansion at level n is a/b
    # analysis shows that at level n+1, we have:
    # (2*b+a) / (a+b)
    # Also, gcd(2*b+a, a+b) = gcd(b, a+b) = gcd(a, b)
    # Thus, gcd is conserved : if we start with a reduced
    # fraction, all fractions will be reduced
    nb = 0
    a = b = 1
    for i in range(nb_exp + 1):
        assert gcd(a, b) == 1
        if len(str(a)) > len(str(b)):
            nb += 1
        a, b = 2 * b + a, a + b
    return nb


def euler58(ratio=0.1):
    """Solution for problem 58."""
    # First analysis in euler28
    # Corners are s*s-3s+3, s*s-2s+2, s*s-s+1, s*s
    # s*s is likely not to be prime
    nb_prime = 0
    for s in itertools.count(3, 2):
        for i in range(1, 4):
            if is_prime(s * s - i * s + i):
                nb_prime += 1
        if nb_prime < ratio * (2 * s - 1):
            return s


def euler59_():
    """Solution for problem 59."""
    pass


def euler60_():
    """Solution for problem 60."""
    # awful C solution on PE https://projecteuler.net/thread=60;page=3#20807
    pass


def euler61_():
    """Solution for problem 61."""
    pass


def euler62(nb_perm=5):
    """Solution for problem 62."""
    cube_perm, l = {}, None
    for i in itertools.count():
        c = i * i * i
        new_l = len(str(c))
        if l != new_l:
            cand = [numbers[0]
                    for c, numbers in cube_perm.items() if len(numbers) == nb_perm]
            if cand:
                return min(cand) ** 3
            cube_perm, l = {}, new_l
        cube_perm.setdefault(sorted_number(c), []).append(i)


def euler63():
    """Solution for problem 63."""
    # 10^(n-1) <= x^n < 10^n
    # (n-1) * log(10) <= n * log(x) < n * log(10)
    # (n-1) * log(10) / n <= log(x) < log(10)
    # exp((n-1) * log(10) / n) <= x < 10
    # and LHS becomes bigger than 9 at n = 22
    return sum(
        10 - ceil(math.exp((n - 1) * math.log(10) / n))
        for n in range(1, 22))


def euler64_():
    """Solution for problem 64."""
    pass


def euler65_():
    """Solution for problem 65."""
    pass


def euler66_():
    """Solution for problem 66."""
    # http://en.wikipedia.org/wiki/Pell's_equation
    pass


def euler67_():
    """Solution for problem 67."""
    pass


def euler68_():
    """Solution for problem 68."""
    # awful C solution on PE https://projecteuler.net/thread=68;page=3#21386
    pass


def euler69(lim=1000000):
    """Solution for problem 69."""
    return max((float(i) / t, i) for i, t in enumerate(totient(lim)) if i)[1]


def euler70(lim=10000000):
    """Solution for problem 70."""
    n, val = lim, lim
    for i, t in enumerate(totient(lim)):
        if i > 1:
            new_val = float(i) / t
            if new_val < val and sorted_number(i) == sorted_number(t):
                n, val = i, new_val
    return n


def euler71_():
    """Solution for problem 71."""
    pass


def euler72(lim=1000000):
    """Solution for problem 72."""
    return sum(totient(lim))


def euler73_():
    """Solution for problem 73."""
    pass


def euler74_():
    """Solution for problem 74."""
    pass


def euler75(lim=1500000):
    """Solution for problem 75."""
    # http://en.wikipedia.org/wiki/Pythagorean_triple
    # a = (m*m-n*n)*d
    # b = 2*m*n*d
    # c = (m*m+n*n)*d
    # p = 2*m*(m+n)*d
    # with : pgcd(m,n)=1 and m < n < 2*m.
    diff_c = [0] * (1 + lim)  # 0: no solution / -1: more than 1 solution
    for m in range(1 + int(math.sqrt(lim / 2))):
        m2 = m * m
        for k in range(m + 1, 2 * m):  # k = n + m
            if gcd(m, k) == 1:
                n = k - m
                n2 = n * n
                mk2 = 2 * m * k
                for d in range(1, 1 + lim // mk2):
                    p = mk2 * d
                    c = (m2 + n2) * d
                    old_res = diff_c[p]
                    if old_res == 0:
                        diff_c[p] = c
                    elif old_res != c:
                        diff_c[p] = -1
    return sum(1 for d in diff_c if d > 0)


def euler76_():
    """Solution for problem 76."""
    pass


def euler77_():
    """Solution for problem 77."""
    pass


def euler78_():
    """Solution for problem 78."""
    pass


def euler79_():
    """Solution for problem 79."""
    pass


def euler80_():
    """Solution for problem 80."""
    pass


def euler81_():
    """Solution for problem 81."""
    pass


def euler82_():
    """Solution for problem 82."""
    pass


def euler83_():
    """Solution for problem 83."""
    pass


def euler84_():
    """Solution for problem 84."""
    pass


def euler85(target=2000000):
    """Solution for problem 85."""
    # nb_rect(i, 1) = 1 (of size i) + 2 (of size i-1) + ... + i (of size 1) = i*(i+1)/2
    # nb_rect(i, j) = nb_rect(i, 1) (of height j) + 2 * nb_rect(i, 1) (of height j-1) + ... + j * nb_rect(i, 1) (of height 1)
    #               = nb_rect(i, 1) * nb_rect(j, 1) = i*(i+1)*j*(j+1)/4
    # Solutions to nb_rect(i, j) = t for a given i are, with delta = 1 + 16 t / [i * (i+1)]
    #   j1 = (-1 - sqrt(delta))/2 < 0 and j2 = (-1 + sqrt(delta))/2 > 0
    # Best approximations are ceil and floor of j2.
    sol_dist, sol_area = target, 0
    for i in itertools.count(1):
        tmp_i = i * (i + 1)
        delta = 1 + target * 16 / tmp_i
        j2 = (-1 + math.sqrt(delta)) / 2
        js = set([math.floor(j2), math.ceil(j2)])
        for j in js:
            val = tmp_i * j * (j + 1) / 4
            dist = math.fabs(target - val)
            if dist < sol_dist:
                sol_dist, sol_area = dist, i * j
        if i > max(js):
            return sol_area


def euler86_():
    """Solution for problem 86."""
    pass


def euler87(lim=50000000):
    """Solution for problem 87."""
    # the biggest prime needed p is such that
    # lim >= p**2 + 2**3 + 2**4
    # p <= sqrt(lim - 24)
    primes = list(primes_up_to(int(math.sqrt(lim - 24))))
    sol = set()
    for a in primes:
        sum_a = a ** 4
        if sum_a > lim:
            break
        for b in primes:
            sum_b = sum_a + b ** 3
            if sum_b > lim:
                break
            for c in primes:
                sum_c = sum_b + c ** 2
                if sum_c > lim:
                    break
                sol.add(sum_c)
    return len(sol)


def euler88_():
    """Solution for problem 88."""
    pass


def euler89_():
    """Solution for problem 89."""
    pass


def euler90():
    """Solution for problem 90."""
    base = 10
    nb_face = 6

    def substitute(n):
        assert 0 <= n < base
        return 6 if n == 9 else n
    # The values we want to display with 9 replaced by 6
    squares = [tuple(substitute(n) for n in divmod(i * i, base)) for i in range(1, base)]
    # Possible dice arrangements with 9 replaced by 6 to make the logic easier
    permutations = (set(substitute(n) for n in p)
                    for p in itertools.combinations(range(base), nb_face))
    # Filtering arrangements as an optimisation
    candidates = (p for p in permutations
                  if all((a in p or b in p) for a, b in squares))
    # Counting the matching arrangements
    return sum(1 for (c1, c2) in itertools.combinations(candidates, 2)
               if all(((a in c1 and b in c2) or (b in c1 and a in c2)) for a, b in squares))


def euler91_bruteforce(size=50):
    """Solution for problem 91.
    Bruteforce is fast enough to get a solution."""
    i = 0
    for p1, q1, p2, q2 in itertools.product(range(size + 1), repeat=4):
        if (p1, p2) < (q1, q2):  # not to count them twice
            low, med, high = sorted(x ** 2 + y ** 2 for x, y in [(p1, p2), (q1, q2), (p1 - q1, p2 - q2)])
            if low and low + med == high:
                i += 1
    return i


def euler91(size=50):
    """Solution for problem 91."""
    # We have 2 kinds of right angles :
    #  - angle in 0 (hypothenus is PQ) : size*size of them
    #  - without any loss of generality, angle in P (p1, p2):
    #       OP is the hypothenus and with H (h1, h2) middle of OP, we have :
    #       |OP| = 2*|OH| = 2*|HP| = 2*|HQ|
    #       p1^2 + p2^2 = 4 * ((q1-h1)^2 + (q2-h2)^2)
    #                   = (2*q1 - p1)^2 + (2*q2 - p2)^2
    #       - If Q is on axis : 2*size*size obvious solutions
    #       - If Q is not on axis, we solve for all p1,p2,q1.
    #   A few optimisations :
    #       - we can assume p1 >= p2 and count solution twice if relevant
    #       - when looking for q1, we can stop early (h1 + |OH|)
    i = 3 * size * size
    for p1, p2 in itertools.combinations_with_replacement(range(size + 1), 2):
        dp = p1 ** 2 + p2 ** 2
        if dp:
            for q1 in range(1, min(size, int((p1 + math.sqrt(dp)) / 2)) + 1):
                square = dp - (2 * q1 - p1) ** 2
                if square >= 0:
                    root = int(math.sqrt(square))
                    if root * root == square:
                        assert p2 % 2 == root % 2
                        for s in set([-root, +root]):
                            q2 = (p2 + s) // 2
                            if 1 <= q2 <= size and (q1, q2) != (p1, p2):
                                i += 1 if p1 == p2 else 2
    return i


def euler92_():
    """Solution for problem 92."""
    pass


def euler93_():
    """Solution for problem 93."""
    pass


def euler94_():
    """Solution for problem 94."""
    pass


def euler95_():
    """Solution for problem 95."""
    pass


def euler96_():
    """Solution for problem 96."""
    pass


def euler97():
    """Solution for problem 97."""
    mod = 10 ** 10
    return (1 + 28433 * pow(2, 7830457, mod)) % mod


def euler98_():
    """Solution for problem 98."""
    pass


def euler99_():
    """Solution for problem 99."""
    pass


def euler100(lim=1000000000000):
    """Solution for problem 100."""
    # P(BB) = (b/t) * ((b-1)/(t-1))
    # P(BB) = 1/2
    # => 2 * b  * (b - 1) = t * (t - 1)
    # https://oeis.org/A046090
    b0, b1 = 1, 3
    r0, r1 = 0, 1
    while True:
        if b0 + r0 > lim:
            return b0
        b0, b1 = b1, 6 * b1 - b0 - 2
        r0, r1 = r1, 6 * r1 - r0


def euler101_():
    """Solution for problem 101."""
    # http://en.wikipedia.org/wiki/Lagrange_polynomial
    # Probably worth creating a polynom module to evaluate/print/interp,etc
    pass


def euler102_():
    """Solution for problem 102."""
    pass


def euler103_():
    """Solution for problem 103."""
    pass


def euler104_(first=True, last=True):
    """Solution for problem 104."""
    digits = sorted('123456789')
    # TODO not efficient enough
    for i, f in enumerate(fibo()):
        s = str(f)
        if (not first or sorted(s[:9]) == digits) \
                and (not last or sorted(s[-9:]) == digits):
            return 1 + i


def euler105_():
    """Solution for problem 105."""
    pass


def euler106_():
    """Solution for problem 106."""
    pass


def euler107_():
    """Solution for problem 107."""
    # http://en.wikipedia.org/wiki/Minimum_spanning_tree
    pass


def increasing_number(n):
    l = str(n)
    return all(l[i] <= l[i + 1] for i in range(len(l) - 1))


def decreasing_number(n):
    l = str(n)
    return all(l[i] >= l[i + 1] for i in range(len(l) - 1))


def bouncy_number(n):
    return not increasing_number(n) and not decreasing_number(n)


def euler108_():
    """Solution for problem 108."""
    pass


def euler109_():
    """Solution for problem 109."""
    pass


def euler110_():
    """Solution for problem 110."""
    pass


def euler111_(nb_dig=4):
    """Solution for problem 111."""
    # Works fine for 4 digits but does not scale well
    lim = 10 ** nb_dig
    low = lim // 10
    prime = sieve(lim)
    m, s = {}, {}
    for i in range(low, lim):
        if prime[i]:
            for dig, nb in collections.Counter(str(i)).items():
                m_val = m.get(dig, 0)
                if nb > m_val:
                    m[dig], s[dig] = nb, [i]
                elif nb == m_val:
                    s[dig].append(i)
    for d in string.digits:
        print(d, m[d], len(s[d]), sum(s[d]))
    print(sum(sum(s[d]) for d in string.digits))


def euler112_(perc=99):
    """Solution for problem 112."""
    nb_bouncy = 0
    for i in itertools.count(1):
        nb_bouncy += 1 if bouncy_number(i) else 0
        if nb_bouncy * 100 == i * perc:
            return i


def euler113(ten_pow=100):
    """Solution for problem 113."""
    # Dynamic programming
    # array with last digit and nb digit
    # Not the cleanest not the most optimal
    nb_incr_number = [1 for i in range(10)]
    nb_decr_number = [0] + [1 for i in range(9)]

    for i in range(ten_pow):
        nb_incr_number = [sum(nb_incr_number[0:i + 1]) for i in range(10)]
        nb_decr_number = [(1 if i else 0) + sum(nb_decr_number[i:10]) for i in range(10)]
    return nb_incr_number[-1] + nb_decr_number[0] - 9 * ten_pow - 1


def euler114(l=50, size_min=3):
    """Solution for problem 114."""
    no_block, block = 1, 0
    nb_comb = []
    for i in range(l):
        # Dynamic programming
        #  Combinations of length n ending with no block are generated taking all combinations of length n-1.
        #  Combinations of length n ending with a block are generated by either :
        #  - adding a block to combinations of length n-1 ending with a block
        #  - starting a new block out of combinations of length n-size_min not ending with a block
        nb_comb.append(no_block)
        no_block, block = (no_block + block), (block + (nb_comb[i - size_min + 1] if i + 1 >= size_min else 0))
    return no_block + block


def euler115(size_min=50, lim=1000000):
    """Solution for problem 115."""
    no_block, block = 1, 0
    nb_comb = []
    for i in itertools.count():
        # Dynamic programming
        #  Combinations of length n ending with no block are generated taking all combinations of length n-1.
        #  Combinations of length n ending with a block are generated by either :
        #  - adding a block to combinations of length n-1 ending with a block
        #  - starting a new block out of combinations of length n-size_min not ending with a block
        nb_comb.append(no_block)
        no_block, block = (no_block + block), (block + (nb_comb[i - size_min + 1] if i + 1 >= size_min else 0))
        if no_block > lim:
            return i


def euler116_():
    """Solution for problem 116."""
    pass


def euler117_():
    """Solution for problem 117."""
    pass


def count_increasing_prime_set_with_digits(digits, start):
    len_dig = len(digits)
    len_str = len(str(start))
    last_dig_for_primes = set('1379')
    for len_perm in range(len_str, 1 + len_dig // 2):
        for comb in itertools.combinations(digits, len_perm):
            rem_dig = digits - set(comb)
            if len_perm == 1 or (sum(int(d) for d in comb) % 3 and
                                 any(d in last_dig_for_primes for d in comb) and
                                 any(d in last_dig_for_primes for d in rem_dig)):
                big_rem = int(''.join(sorted(rem_dig, reverse=True)))
                cand = [n for n in (int(''.join(p)) for p in itertools.permutations(comb)) if start < n < big_rem and is_prime(n)]
                if cand:
                    mini = min(cand)
                    lst = list(count_increasing_prime_set_with_digits(rem_dig, mini))
                    for c in cand:
                        nb = sum(num for small, num in lst if c < small)
                        if nb:
                            yield (c, nb)

    if len_dig == 1 or (sum(int(d) for d in digits) % 3 and any(d in last_dig_for_primes for d in digits)):
        for c in [n for n in (int(''.join(p)) for p in itertools.permutations(digits)) if start < n and is_prime(n)]:
            yield (c, 1)


def euler118():
    """Solution for problem 118."""
    return sum(s for _, s in count_increasing_prime_set_with_digits(set('123456789'), 0))


def yield_interesting_119(base=10):
    """Generate interesting numbers for problem 119."""
    for nb_digits in itertools.count(2):
        min_val = base ** (nb_digits - 1)
        max_val = min_val * 10
        candidates = []
        for dig_sum in range(2, (base - 1) * nb_digits):
            n = dig_sum
            while n < min_val:
                n *= dig_sum
            while n < max_val:
                if sum(int(d) for d in str(n)) == dig_sum:
                    candidates.append(n)
                n *= dig_sum
        for c in sorted(candidates):
            yield c


def euler119(n=30):
    """Solution for problem 119."""
    return nth(yield_interesting_119(), n - 1)


def euler120(beg=3, end=1000):
    """Solution for problem 120."""
    # A(n) = (a-1)^n + (a+1)^n
    #      = sum(k=0..n, (n, k) a^(k) (-1)^(n-k)) + sum(k=0..n, (n, k) a^(k) (1)^(n-k))
    #      = sum(k=0..n, (n, k) a^(k) (-1)^(n-k)) + sum(k=0..n, (n, k) a^(k) (1)^(n-k))
    #      = sum(k=0..n, (n, k) a^(k) [ (-1)^(n-k) + (1)^(n-k)] )
    #      = sum(k=0..n if (n-k) is even, 2 (n, k) a^(k))
    # B(n) = A(n)                                            [a^2]
    #      = sum(k=0..n if (n-k) is even, 2 (n, k) a^(k))    [a^2]
    # Modulo a^2, all terms with k >= 2 are O
    # B(n) = (2 (n, 0) a^0 if n is even) + (2 (n, 1) a^1 if n-1 is even) [a^2]
    # B(n) = (2 if n is even else 2*n*a)                                 [a^2]
    # C(n) = 2*n*a [a^2]
    #    C(n+ia) = 2*(n+ia)*a      [a^2]
    #            = 2*n*a + 2*i*a^2 [a^2]
    #            = 2*n*a           [a^2]
    #            = C(n)            [a^2]
    #       => Period a
    # If a is even (and > 2):
    #    Period is actually a/2
    #    0    < 2*a  < 4*a  < ... < a^2 - 2a   > 0
    #    C(0) < C(1) < C(2) < ... < C(a/2 - 1) > C(a/2) = C(0) = 0
    #    so the max for C on |N is a^2 - 2a reached for n = a/2 - 1 (+ka/2)
    #    so the max for B on odd numbers is a^2 - 2a
    #       (adding a/2 to n to reach an odd number if needed)
    #    max(B) = max(2 % a^2, a^2 - 2a)
    # If a is odd (and > 1):
    #    0    < 2*a  < 4*a  < ... < a^2 - a    > a
    #    C(0) < C(1) < C(2) < ... < C((a-1)/2) > C((a+1)/2)
    #    so the max for C is reached for n = (a-1)/2 on the first half-period
    #    a          < 3*a  < 5*a  < a^2 - 2a > 0
    #    C((a+1)/2) < ...         < C(a-1)   > C(a) = C(0) = 0
    #    so the max for C is reached for n = (a-1) on the second half-period
    #    the max for C on |N is a^2 - a reached for n = (a-1)/2 (+ka)
    #    so the max for B on odd numbers is a^2 - a
    #       (adding a to n to reach an odd number if needed)
    #    max(B) = max(2 % a^2, a^2 - a)
    #
    # Finally for all a bigger than 2, max(B) = a^2 - a if a if even else a^2 - 2a
    # These could be added with a proper summation formula but it is not worth
    # the pain.
    assert beg > 2
    return sum(a * (a - (1 if a % 2 else 2)) for a in range(beg, end + 1))


def euler121_():
    """Solution for problem 121."""
    pass


def euler122_():
    """Solution for problem 122."""
    pass


def euler123_():
    """Solution for problem 123."""
    # C solution on PE https://projecteuler.net/thread=123;page=3#21789
    pass


def euler124(lim=100000, n=10000):
    """Solution for problem 124."""
    rad = [mult(div) for div in prime_divisors_sieve(lim)]
    rev_rad = {}
    for i, r in enumerate(rad):
        rev_rad.setdefault(r, []).append(i)
    return nth((i for rad in sorted(rev_rad.keys()) for i in rev_rad[rad]), n)


def euler125(lim=100000000):
    """Solution for problem 125."""
    # Pointless analysis I keep for later:
    # S(n) = 1^2 + 2^2 + ... + n^2 = n * (n+1) * (2n+1) / 6
    # S(k, n) = k^2 + (k+1)^2 + ... + n^2
    #         = 1^2 + 2^2 + ... + n^2 - [1^2 + 2^2 + ... (k-1)^2]
    #         = S(n) - S(k-1)
    # S(k, k+n) = S(k+n) - S(k-1)
    #           = [(k+n) * (k+n+1) * (2k+2n + 1) - (k-1) * k * (2k-1)] / 6
    #           = (6 k^2 n + 6 k^2 + 6 k n^2 + 6kn + 2 n^3 + 3 n^2 + n) / 6
    #           = k^2 n + k^2 + k n^2 + kn + n^3 / 3 + n^2 / 2 + n/6
    # Upper bound for k and n - could be improved with analysis above
    # but this is fast enough
    sqrt_lim = int(math.sqrt(lim))
    found = set()
    for k in range(1, sqrt_lim + 1):
        s = k * k
        for n in range(k + 1, sqrt_lim + 1):
            s += n * n
            if s > lim:
                break
            else:
                st = str(s)
                if st == st[::-1]:
                    found.add(s)
    return sum(found)


def euler126_():
    """Solution for problem 126."""
    pass


def euler127(lim=120000):
    """Solution for problem 127."""
    # 1) If a, b, c are relatively prime, rad(abc) = rad(a) * rad(b) * rad(c)
    # 2) gcd(a, b) == 1 <=> gcd(rad(a), rad(b)) == 1
    # 3) 1 = gcd(a, c) = gcd(a, c - a) = gcd(a, b) = gcd(a+b, b) = gcd(c, b)
    # 4) As c > rad(abc) > rad(c), we must have rad(c) < c which implies that
    #       2 * rad(c) <= c
    # Technique is the following :
    #  - create a cache simulating the rad function
    #  - create a map from rad values to the (ordered) list of values
    #  - iterate over increasing rad(c) (we can stop early)
    #    - iterate over decreasing c (we can stop early)
    #      - iterate over increasing rad(a) (we can stop early)
    #        - iterate over increasing a (we can stop early)
    #          - check a, compute b, check b
    rad = [mult(div) for div in prime_divisors_sieve(lim - 1)]
    rev_rad = {}
    for i, r in enumerate(rad):
        rev_rad.setdefault(r, []).append(i)
    sorted_rad = sorted(rev_rad.keys())
    s = 0
    for radc in sorted_rad:
        if 2 * radc > lim:
            break
        for c in reversed(rev_rad[radc]):
            if 2 * radc >= c:
                break
            for rada in sorted_rad:
                if 2 * rada * radc >= c:
                    break
                if gcd(rada, radc) == 1:
                    for a in rev_rad[rada]:
                        b = c - a
                        if a >= b:
                            break
                        if b < c and rada * rad[b] * radc < c:
                            s += c
    return s


def euler128_():
    """Solution for problem 128."""
    pass


def euler129_():
    """Solution for problem 129."""
    pass


def euler130_():
    """Solution for problem 130."""
    pass


def euler131(lim=1000000):
    """Solution for problem 131."""
    # m^3 = n^3 + p*n^2 = n*n*(n+p)
    # 1) Proof that this is equivalent to n and (n+p) being perfect cube
    # Let m     = prod(pi^ei) be the prime decomposition of m
    #     n     = prod(pi^fi)                               n
    #     (n+p) = prod(pi^gi)                               (n+p)
    # Let's prove that n and (n+p) must be perfect cube ie gi and fi
    # are all multiple of 3.
    # m^3 = prod(pi^(3*ei))      is the prime decomposition of m^3
    #     = prod(pi^(2*fi + gi))
    # For all primes pi, 3*ei = 2*fi + gi
    # For any "real prime factor" (ie, ei != 0):
    #  - if fi = 0 => gi is a multiple of 3 (and so is fi)
    #  - if gi = 0 => fi is a multiple of 3 (and so is gi)
    #  - if fi != 0 and gi != 0:
    #       => pi divides n and n+p
    #       => pi divides p (and p is prime)
    #       => pi = p
    #       => n = k*p
    #       => m^3 = k*p*k*p*(k*p+1*p) = p^3 * (k*k*(k+1))
    #       Dividing by p^3 on both side, we'd need k*k*(k+1) to be a
    #       perfect cube but its obviously between two consecutives cubes :
    #       k^3 and (k+1)^3 so this is not possible.
    # Therefore, for each prime divisor pi, pi has a power divisible by 3 in
    # n and n+p so n and n+p are both cubes.
    # In the other direction, if n and (n+p) are perfect cube, n*n*(n+p)
    # obviously is.
    # 2) Proof that p = 3*b^2 + 3*b + 1 with b positive integer
    # If n = b^3 and n+p = a^3 with a > b > 0
    # Then, p = n+p - n = a^3 - b^3 = (a-b) (a^2 + ab + b^2)
    # Because p is prime, we must have : a - b = 1 <=> a = b + 1
    # So p = (b+1)^2 + (b+1)*b + b^2 = b^2 + 2b + 1 + b^2 + b + b^2
    #      = 3*b^2 + 3*b + 1
    # Any *prime* value 3*b^2 + 3*b + 1 would be a solution
    # (with n = b^3 and m = b*b*a = b*b*(b+1)).
    count = 0
    for b in itertools.count():
        p = 3 * b * b + 3 * b + 1
        if p > lim:
            break
        if is_prime(p):
            count += 1
    return count


def euler132_():
    """Solution for problem 132."""
    pass


def euler133_():
    """Solution for problem 133."""
    pass


def euler134_():
    """Solution for problem 134."""
    pass


def euler135_():
    """Solution for problem 135."""
    pass


def euler136_():
    """Solution for problem 136."""
    pass


def euler137_():
    """Solution for problem 137."""
    pass


def euler138_():
    """Solution for problem 138."""
    pass


def euler139_():
    """Solution for problem 139."""
    pass


def euler140_():
    """Solution for problem 140."""
    pass


def euler141(lim=10 ** 12):
    """Solution for problem 141."""
    # m^2 = n = d*q + r with 0 <= r < d
    # Values can ke written c, c*k, c*k*k with k > 0 (values are all positive)
    # We can assume k > 1 (or we take 1/k)
    # We have different possible orders :
    # d, q, r and q, d, r are not possible because r < d
    # d, r, q and q, r, d are not possible because it leads to r*r = q*d which means n = r*(r+1), not a square
    # r, q, d and r, d, q are somewhat symetric cases as q and d have similar roles
    # Consecutives terms will be : r, rk, rkk
    # k is rational because a and a*k are both integers
    # k can be expressed in an irreductible form k = a/b with a > b > 0 and gcd(a, b) == 1
    # Consecutives terms will be : r, r*a/b, r*a^2/b^2
    # Because all terms (especially the last one) are integers, we must have b^2 dividing r
    # So r = c*b*b and other terms are a*b*c and a*a*c
    # And we have :
    # m^2 = n = dq + r = a^3*b*c^2 + b^2*c
    sol = set()
    for a in itertools.count(2):
        a3 = a * a * a
        if a3 >= lim:
            break
        for b in range(1, a):
            b2 = b * b
            if a3 * b + b2 >= lim:
                break
            if gcd(a, b) == 1:
                for c in itertools.count(1):
                    n = a3 * b * c * c + b2 * c
                    if n >= lim:
                        break
                    sqrt = int(math.sqrt(n))
                    if sqrt * sqrt == n:
                        sol.add(n)
    return sum(sol)


def euler142_():
    """Solution for problem 142."""
    pass


def euler143_():
    """Solution for problem 143."""
    pass


def euler144_():
    """Solution for problem 144."""
    pass


def euler145_():
    """Solution for problem 145."""
    pass


def euler146_():
    """Solution for problem 146."""
    pass


def euler147_():
    """Solution for problem 147."""
    pass


def euler148_():
    """Solution for problem 148."""
    pass


def euler149_():
    """Solution for problem 149."""
    pass


def euler150_():
    """Solution for problem 150."""
    pass


def euler151_():
    """Solution for problem 151."""
    pass


def euler152_():
    """Solution for problem 152."""
    pass


def euler153_():
    """Solution for problem 153."""
    pass


def euler154_():
    """Solution for problem 154."""
    pass


def euler155_():
    """Solution for problem 155."""
    pass


def euler156_():
    """Solution for problem 156."""
    pass


def euler157_():
    """Solution for problem 157."""
    pass


def euler158_():
    """Solution for problem 158."""
    pass


def euler159_():
    """Solution for problem 159."""
    pass


def euler160_():
    """Solution for problem 160."""
    pass


def euler161_():
    """Solution for problem 161."""
    pass


def euler162_():
    """Solution for problem 162."""
    pass


def euler163_():
    """Solution for problem 163."""
    pass


def euler164():
    """Solution for problem 164."""
    # Dynamic programming
    nb = {i: 1 for i in range(1, 10)}
    for _ in range(2, 21):
        new = collections.defaultdict(int)
        for i, n in nb.items():
            i1, i2 = divmod(i, 10)
            for i3 in range(10 - i1 - i2):
                new[10 * i2 + i3] += n
        nb = new
    return sum(nb.values())


def euler165_():
    """Solution for problem 165."""
    # http://en.wikipedia.org/wiki/Line_segment_intersection
    pass


def euler166_():
    """Solution for problem 166."""
    pass


def euler167_():
    """Solution for problem 167."""
    pass


def euler168_():
    """Solution for problem 168."""
    pass


def euler169_():
    """Solution for problem 169."""
    pass


def euler170_():
    """Solution for problem 170."""
    pass


def euler171_():
    """Solution for problem 171."""
    pass


def euler172_():
    """Solution for problem 172."""
    pass


def euler173_():
    """Solution for problem 173."""
    # decent looking C on PE https://projecteuler.net/thread=173;page=2#21231
    pass


def euler174_():
    """Solution for problem 174."""
    pass


def euler175_():
    """Solution for problem 175."""
    pass


def euler176_():
    """Solution for problem 176."""
    pass


def euler177_():
    """Solution for problem 177."""
    pass


def euler178_():
    """Solution for problem 178."""
    pass


def euler179_():
    """Solution for problem 179."""
    # braceless C on PE https://projecteuler.net/thread=179;page=4#21050
    pass


def euler180_():
    """Solution for problem 180."""
    pass


def euler181_():
    """Solution for problem 181."""
    pass


def euler182_():
    """Solution for problem 182."""
    pass


def k_max_for_euler183(n):
    # P(x) = (N/x)^x has a maximum in x = N/e
    m = n / math.e
    return max(math.floor(m), math.ceil(m), key=lambda val: val * math.log(n / val))  # keeps ordering


def euler183(lim=10000):
    """Solution for problem 183."""
    def fraction_is_finite(a, b, base):
        # computing powers does not change "finiteness"
        # (a/b)^n is as "finite" as (a/b)
        # which is as "finite" as :
        # (a/gcd(a,b)) / (b/gcd(a,b)) or
        # 1 / (b/gcd(a,b))
        b /= gcd(a, b)
        for d in range(2, base):
            if base % d == 0:
                while b % d == 0:
                    b /= d
        return b == 1
    return sum(n * (-1 if fraction_is_finite(n, k_max_for_euler183(n), 10) else 1)
               for n in range(5, lim + 1))


def euler184_():
    """Solution for problem 184."""
    pass


def euler185_():
    """Solution for problem 185."""
    pass


def euler186_():
    """Solution for problem 186."""
    pass


def euler187_():
    """Solution for problem 187."""
    # slow C solution on PE https://projecteuler.net/thread=187#20276
    pass


def euler188_():
    """Solution for problem 188."""
    pass


def euler189_():
    """Solution for problem 189."""
    pass


def euler190_():
    """Solution for problem 190."""
    # C solution on PE https://projecteuler.net/thread=190#21099
    pass


def euler191(nb_days=4):
    """Solution for problem 191."""
    # 0 late : a, b, c for 0, 1, 2 consecutive absences
    # 1 late : d, e, f for 0, 1, 2 consecutive absences
    a, b, c, d, e, f = 1, 0, 0, 0, 0, 0
    for _ in range(nb_days + 1):  # 1 more iteration to have the res in d
        a, b, c, d, e, f = a + b + c, a, b, a + b + c + d + e + f, d, e
    return d


def euler192_():
    """Solution for problem 192."""
    pass


def euler193_():
    """Solution for problem 193."""
    pass


def euler194_():
    """Solution for problem 194."""
    pass


def euler195_():
    """Solution for problem 195."""
    pass


def euler196_():
    """Solution for problem 196."""
    pass


def euler197_():
    """Solution for problem 197."""
    pass


def euler198_():
    """Solution for problem 198."""
    pass


def euler199_():
    """Solution for problem 199."""
    pass


def euler200_():
    """Solution for problem 200."""
    pass


def euler201_():
    """Solution for problem 201."""
    pass


def euler202_():
    """Solution for problem 202."""
    pass


def euler203_():
    """Solution for problem 203."""
    pass


def euler204_():
    """Solution for problem 204."""
    pass


def euler205_():
    """Solution for problem 205."""
    pass


def euler206_():
    """Solution for problem 206."""
    pass


def euler207_():
    """Solution for problem 207."""
    pass


def euler208_():
    """Solution for problem 208."""
    pass


def euler209_():
    """Solution for problem 209."""
    pass


def euler210_():
    """Solution for problem 210."""
    pass


def euler211_():
    """Solution for problem 211."""
    pass


def euler212_():
    """Solution for problem 212."""
    pass


def euler213_():
    """Solution for problem 213."""
    pass


def euler214(lim=40000000, length=25):
    """Solution for problem 214."""
    tot = totient(lim)
    chains = [0] * (lim + 1)
    for i, t in enumerate(tot):
        chains[i] = 1 + (chains[t] if i > 1 else 0)
    return sum(i for i, (l, t) in enumerate(zip(chains, tot)) if l == length and i == t + 1)


def euler215(L=32, H=10):
    """Solution for problem 215."""
    # I have no idea how this was supposed to work
    # but there must be a better/faster way
    crack = {0: [[]]}
    for l in range(1, L + 1):
        crack[l] = [i + ([l - s] if l - s > 0 else []) for s in [2, 3] for i in crack.get(l - s, [])]
    rows = crack[L]
    compat = [(i, j) for (i, r), (j, s) in itertools.combinations(enumerate(rows), 2) if all(e not in r for e in s)]

    comb = {i: 1 for i in range(len(rows))}
    for h in range(1, H):
        new_comb = collections.defaultdict(int)
        for i, j in compat:
            new_comb[j] += comb[i]
            new_comb[i] += comb[j]
        comb = new_comb
    return sum(comb.values())


def euler216_():
    """Solution for problem 216."""
    pass


def euler217_():
    """Solution for problem 217."""
    pass


def euler218_():
    """Solution for problem 218."""
    pass


def euler219_():
    """Solution for problem 219."""
    pass


def euler220_():
    """Solution for problem 220."""
    pass


def euler221_():
    """Solution for problem 221."""
    pass


def euler222_():
    """Solution for problem 222."""
    pass


def euler223_():
    """Solution for problem 223."""
    pass


def euler224_():
    """Solution for problem 224."""
    pass


def euler225(num=124):
    """Solution for problem 225."""
    def has_mod_property(n):
        s = set()
        t1 = t2 = t3 = 1
        while True:
            if t3 == 0:
                return False
            if ((t1, t2, t3)) in s:
                return True
            s.add((t1, t2, t3))
            t1, t2, t3 = t2, t3, (t1 + t2 + t3) % n

    return nth((i for i in itertools.count(1, 2) if has_mod_property(i)), num - 1)


def euler226_():
    """Solution for problem 226."""
    pass


def euler227_():
    """Solution for problem 227."""
    pass


def euler228_():
    """Solution for problem 228."""
    pass


def euler229_():
    """Solution for problem 229."""
    pass


def euler230_():
    """Solution for problem 230."""
    pass


def euler231_():
    """Solution for problem 231."""
    pass


def euler232_():
    """Solution for problem 232."""
    pass


def euler233_():
    """Solution for problem 233."""
    pass


def euler234_():
    """Solution for problem 234."""
    pass


def euler235_(x=900, y=-3, n=5000, v=-600000000000):
    """Solution for problem 235."""
    # http://en.wikipedia.org/wiki/Arithmetico-geometric_sequence
    pass


def euler236_():
    """Solution for problem 236."""
    pass


def euler237_():
    """Solution for problem 237."""
    pass


def euler238_():
    """Solution for problem 238."""
    pass


def euler239_():
    """Solution for problem 239."""
    pass


def euler240_():
    """Solution for problem 240."""
    pass


def euler241_():
    """Solution for problem 241."""
    pass


def euler242_():
    """Solution for problem 242."""
    pass


def euler243(num=15499, den=94744):
    """Solution for problem 243."""
    # R(d) = phi(d) / (d-1)
    # Small values can be obtained by multiplying primes : 2*3*5*...
    # From there a bit of guessing, hard to write automated code for this.
    return 892371480


def euler301(max_pow_of_2=30):
    # Computing whether a Nim game can be won or not (and thus, computing X),
    # corresponds to computing the xor of the different heaps size (aka nim-sum)
    # See https://en.wikipedia.org/wiki/Nim for more details
    # In particular, X(n1, n2, n3) = 0 iff n1 ^ n2 ^ n3 = 0
    # For instance, X(1, 2, 3) = 0
    #
    # This can also be seen as: X(n1, n2, n3) = 0 iff n1 ^ n2 = n3
    # In the case of the problem with X(n, 2*n, 3*n), we want to find values where:
    # n ^ (2*n) = (3*n) = n + (2*n).
    # This corresponds to finding the cases (n, 2*n) where xor (with no carries)
    # acts just like the usual addition (with carries): for a given bit position,
    # the value can't be 1 in both n and 2n.
    # Because the binary representation for (2*n) is the same as the one for n but
    # shifted by a bit, X(n, 2*n, 3*n) = 0 iff there is no consecutive 1 in the
    # binary representation of n.
    # We can count such occurences of a binary representation of length l:
    #  nb(l=0) = 1 ("")
    #  nb(l=1) = 2 ("0", "1")
    #  nb(l=2) = 3 ("00", "01", "10")
    #  nb(l+2) = nb(l+1)    <- add a '0' to a smaller string
    #            + nb(l)    <- add a '01' to a smaller string
    #  => nb(l) = Fibo(l)
    return nth(fibo(1, 2), max_pow_of_2)


def euler387(nb_zeros=14):
    """Solution for problem 387."""
    # TODO: Most of the time is spent checking primeness of big numbers
    # Using a Rabin Miller approach fixes this
    base = 10  # Not so much of a variable due to various optimisations
    new_rt_harshald = range(1, base)
    strong_rt_harshald_primes = []
    for it in range(nb_zeros - 2):
        prev_rt_harshald, new_rt_harshald = new_rt_harshald, list()
        for p in prev_rt_harshald:
            b = base * p
            for i in range(b, b + base):
                q, r = divmod(i, sum_digit(i))
                if r == 0:
                    # assert len(str(i)) == it+2
                    new_rt_harshald.append(i)
                    if is_prime(q):
                        b2 = base * i
                        for j in [b2 + 1, b2 + 3, b2 + 7, b2 + 9]:
                            if is_prime(j):
                                # assert len(str(j)) == it+3
                                strong_rt_harshald_primes.append(j)
    return sum(strong_rt_harshald_primes)


def euler244_():
    """Solution for problem 244."""
    pass


def euler245_():
    """Solution for problem 245."""
    pass


def euler246_():
    """Solution for problem 246."""
    pass


def euler247_():
    """Solution for problem 247."""
    pass


def euler248_():
    """Solution for problem 248."""
    pass


def euler249_():
    """Solution for problem 249."""
    pass


def euler250_():
    """Solution for problem 250."""
    pass


def euler265_():
    """Solution for problem 265."""
    pass  # Solution on bitbucket to be copied


def euler303_():
    """Solution for problem 303."""
    pass  # Solution on bitbucket to be copied


def euler491():
    """Solution for problem 491."""
    # A number is divisible by 11 if and only the difference between the sum
    # of the digits of odd index and the sum of the digits of even index is
    # divisible by 11.
    # A solution is to generate all distinct combinations of 10 digits and
    # after checking that it leads to a difference divisible by 11, compute
    # the corresponding contribution.
    # For each distinct comb c, if comp is the complement of c in double-pand,
    # contribution is :
    # C = number_of_valid_perm(c) * number_of_valid_perm(comp)
    #   = nb_perm_with_rep_without_lead_0(c) * nb_perm_with_rep(comp)
    #     len(c) - nb(0, c)       len(c) !        len(comb) !
    #   = ----------------- * ----------------- * ------------------
    #           len(c)        prod(multipli_c!)   prod(multipl_comp!)
    # because each elements appears at most twice and number of dup in
    #  c is the number of dup in comp
    #      nb(non-0, c) * (len(c) - 1) ! * len(c) !
    # C =  ----------------------------------------
    #        2 ^ (2 * nb_duplicates(c))
    pand = list(range(10))
    s = 2 * sum(pand)
    comb_len = len(pand)
    prod_fact = math.factorial(comb_len) * math.factorial(comb_len - 1)
    total = 0
    # Next 6 lines are just a tedious but efficient way to generate distinct
    # combinations. "for c in itertools.combinations(numbers, comb_len):"
    # would work too but would need to have duplicates removed (with a set for
    # instance). However, it would be pretty tedious to adapt if the format
    # was to ever change.
    # It first selects the distinct digits that are used, then add the
    # duplicates.
    for l in range(comb_len + 1):  # first 5 iterations are pointless but...
        nb_dup = comb_len - l
        pow_two = 2 ** (2 * nb_dup)
        for dist_dig in itertools.combinations(pand, l):
            for dup in itertools.combinations(dist_dig, nb_dup):
                c = dist_dig + dup
                if (2 * sum(c) - s) % 11 == 0:
                    total += len([x for x in c if x]) * prod_fact // pow_two
    return total


def euler500_(pow_of_two=500500, mod=500500507):
    """Solution for problem 500."""
    # nb_divisors(p_i^k_i * ... p_n^k_n) = (1+k_i) * .. (1+k_n)
    # when nb_divisors is a power of 2, all terms must be powers of 2.
    # we are looking for smallest value : p_1^k_1 * ... p_n^k_n
    # such that 1+k_i are powers of 2
    # so powers (k_i) can be written 2^n - 1.
    prime_gen = yield_primes()
    last_prime = next(prime_gen)
    # heap = [(math.log(last_prime), last_prime, 1)]
    heap = [(last_prime, last_prime, 1)]
    for _ in range(pow_of_two):
        val, prime, power = heapq.heappop(heap)
        # heapq.heappush(heap, (2 * power * math.log(prime), prime, 2*power))
        heapq.heappush(heap, (val ** 2, prime, 2 * power))
        if prime == last_prime:
            last_prime = next(prime_gen)
            # heapq.heappush(heap, (math.log(last_prime), last_prime, 1))
            heapq.heappush(heap, (last_prime, last_prime, 1))
    return mult(pow(prime, power - 1, mod) for val, prime, power in heap) % mod


def euler684(lim=90, mod=1000000007, base=10):

    def s_naive(n):
        q, r = divmod(n, base-1)
        return int(str(r) + q*str(base-1))

    assert s_naive(10) == 19

    def s(n):
        q, r = divmod(n, base-1)
        return (r+1) * pow(base, q, mod) - 1

    for i in range(30):
        assert s(i) == s_naive(i)

    def S_naive(k):
        return sum(s(i) for i in range(1, k+1)) % mod

    assert S_naive(20) == 1074

    def S(k):
        # Optimisation to group numbers that would be linked
        # to the same power of base in s(n)
        q, r = divmod(k, base-1)
        base_pow_q = pow(base, q, mod)

        # For the q first powers considered, there are
        # 'base-1' terms that will be associated to the same power of base
        # with different multiplicative coeff
        # Sum (2 + 3 + ... + base) = base * (base+1) / 2 - 1
        sum_coef_1 = base * (base+1) // 2 - 1
        # The sum of these q first powers: base^0 + base^1 + .. + base^(q-1)
        # is: (base^q - 1) / (base - 1)
        inv_base_mod = pow(base - 1, -1, mod)
        sum_pow = inv_base_mod * (base_pow_q - 1) % mod

        # For the biggest power considered, there are 'r' remaining
        # terms that will be associated with different multiplicative coeff
        # Sum(2 + .. + r+1) = (r+1) * (r+2) / 2 - 1
        sum_coef_2 = (r+1) * (r+2) // 2 - 1
        s = sum_coef_1 * sum_pow + sum_coef_2 * base_pow_q - q * (base - 1) - r
        return s % mod

    for i in range(2000):
        assert S(i) == S_naive(i)

    fib = fibo(0, 1)
    next(fib)
    next(fib)
    return sum(S(f) for _, f in zip(range(2, lim+1), fib)) % mod


def euler686(L=123, n=678910, base=10):
    # 2^j starts with L
    # There is a i, such that:
    # A = L*10^i <= 2^j < (L+1)*10^i = B
    # log(A) <= j*log(2) < log(B)
    # log(A)/log(2) <= j < log(B) / log(2)
    log_base = math.log(base)
    log_l = math.log(L)
    log_lp1 = math.log(L+1)
    log_2 = math.log(2)
    lst = []
    for i in itertools.count(0):
        i_log_base = i * log_base
        ceil_lhs = math.ceil((log_l + i_log_base) / log_2)
        ceil_rhs = math.ceil((log_lp1 + i_log_base) / log_2)
        for v in range(ceil_lhs, ceil_rhs):
            lst.append(v)
            if len(lst) == n:
                return v


def can_be_split_in_sum(digits, target, base = 10):
    """Check if target can be reached by summing part of digits."""
    # Examples:
    # 81 9 -> true (because 8+1)
    # 6724 82 -> true (because 6+72+4)
    # 8281 91 -> true (because 8+2+81)
    # 9801 99 -> true (because 98+0+1)
    # 100 1 -> true (because 1+0)
    if digits < target:
        return False
    if digits == target:
        return True
    power = base
    while power <= digits:
        # power increases, tail increases, head decreases
        head, tail = divmod(digits, power)
        assert head
        if target < tail:
            return False
        if can_be_split_in_sum(head, target - tail):
            return True
        power *= base
    assert divmod(digits, power) == (0, digits)
    return False


def euler700(mul=1504170715041707, mod=4503599627370517):
    def naive():
        lst = []
        latest = None
        for i in range(1, mod):
            v = mul * i % mod
            if latest is None or v < latest:
                latest = v
                lst.append(latest)
        return lst
    if mul < 100000:
        print(mul, mod, naive())

    assert gcd(mul, mod) == 1
    inv = modinv(mul, mod)
    print(mul, mod, inv)
    latest = mod
    latest_i = 0
    i = 1
    while True:
         v = i * inv % mod
         assert i == v*mul%mod
         if v < latest:
             print(i, v, i - latest_i, latest-v)
             i=2*i-latest_i
             latest = v
             latest_i = i
             if v == 1:
                 break
         else:
             i+=1


def euler719(n=10**12):
    """Solution for problem 719."""
    return sum(i*i
               for i in range(2, 1 + int(math.sqrt(n)))
               if can_be_split_in_sum(i*i, i))

def euler759():
    """Solution for problem 759."""
    pass

def euler788(power_input=2022, base=10, modulo=1000000007):
    """Solution for problem 759."""

    def bruteforce_analysis():
        """Bruteforce analysis to confirm hypothesis."""
        elements = dict()
        for i in range(1, base**power_input):
            s = str(i)
            N = len(s)
            c = collections.Counter(s)
            d, n = c.most_common()[0]
            if 2*n > N:
                elements.setdefault((N,n), []).append(i)
        return sum(len(v) for v in elements.values())

    # How many many numbers of length N have exactly n digits equal to d (with 2n > N)
    #  - numbers starting with d:
    #     * d can be any value in [1, base-1]
    #     * remaining (n-1) d can be combined in (N-1 choose n-1) ways
    #     * remaining digits can be any other numbers
    #     This leads to: (base-1) * (base-1)**(N-n) * C(N-1, n-1) distinct options
    #
    #  - numbers not starting with d:
    #     * first digit must not be 0
    #     * d can be any value in [0, base-1] except the first digit
    #     * d can be combined in (N-1 choose n) ways
    #     * remaining digits can be any other numbers
    #     This leads to: (base-1) * (base-1) * (base-1)**(N-n-1) * C(N-1, n) distinct options
    #
    #  The result is:
    #  R =  (base-1) * (base-1)**(N-n) * C(N-1, n-1) + (base-1) * (base-1) * (base-1)**(N-n-1) * C(N-1, n)
    #    =  (base-1)**(N-n+1) * C(N-1, n-1) + (base-1)**(N-n+1) * C(N-1, n)
    #    =  (base-1)**(N-n+1) * [ C(N-1, n-1) + * C(N-1, n)]
    #    =  (base-1)**(N-n+1) * C(N, n)
    #
    # TODO: Summing values of R over different values of N and n could be optimised
    # by reorganising terms and use the hockey stick identity (?):
    #
    # S = Sum(N from 1 to input, n from 1+N/2 to N, (base-1)**(N-n+1) * C(N, n))
    #   = Sum(N from 1 to input, rem from 0 to N/2, (base-1)**(rem+1) * C(N, N-rem))
    #   = Sum(N from 1 to input, rem from 0 to N/2, (base-1)**(rem+1) * C(N, rem))
    #   = Sum(N from 1 to input, rem from 0 to N/2, (base-1)**(rem+1) * C(N, rem))
    R = 0
    for rem in range((1+power_input)//2):
        power = pow(base-1, rem+1, modulo)
        for N in range(2*rem+1, 1+power_input):
            R += power * math.comb(N, rem)
    return R % modulo

tests = [
    (euler1, (10,), 23),
    (euler1, (), 233168),
    (euler2, (), 4613732),
    (euler3, (13195,), 29),
    (euler3, (), 6857),
    (euler4, (2,), 9009),
    (euler4, (), 906609),
    (euler5, (10,), 2520),
    (euler5, (), 232792560),
    (euler6, (10,), 2640),
    (euler6, (), 25164150),
    (euler7, (6,), 13),
    (euler7, (), 104743),
    (euler8, (4,), 5832),
    (euler8, (), 23514624000),
    (euler9, (), 31875000),
    (euler10, (10,), 17),
    (euler10, (), 142913828922),
    (euler12, (5,), 28),
    (euler12, (), 76576500),
    (euler14, (), 837799),
    (euler15, (2, 2), 6),
    (euler15, (), 137846528820),
    (euler16, (15,), 26),
    (euler16, (), 1366),
    (euler17, (5, ), 19),
    (euler17, (), 21124),
    (euler19, (), 171),
    (euler20, (10,), 27),
    (euler20, (), 648),
    (euler21, (), 31626),
    (euler22, (), 871198282),
    (euler23, (), 4179871),
    (euler24, (), 2783915460),
    (euler25, (3,), 12),
    (euler25, (), 4782),
    (euler26, (11,), 7),
    (euler26, (), 983),
    (euler27, (), -59231),
    (euler28, (5,), 101),
    (euler28, (), 669171001),
    (euler29, (5, 5), 15),
    (euler29, (), 9183),
    (euler30, (), 443839),
    (euler31, (), 73682),
    (euler32, (), 45228),
    (euler33, (), 100),
    (euler34, (), 40730),
    (euler35, (2,), 13),
    (euler35, (), 55),
    (euler36, (), 872187),
    (euler37, (), 748317),
    (euler38, (), 932718654),
    (euler39, (), 840),
    (euler40, (), 210),
    (euler41, (), 7652413),
    (euler43, (), 16695334890),
    (euler44, (), 5482660),
    (euler45, (), 1533776805),
    (euler46, (), 5777),
    (euler47, (2,), 14),
    (euler47, (3,), 644),
    (euler47, (), 134043),
    (euler48, (10, 10), 405071317),
    (euler48, (), 9110846700),
    (euler49, (1,), None),
    (euler49, (2,), None),
    (euler49, (3,), None),
    (euler49, (), 296962999629),
    (euler50, (100,), 41),
    (euler50, (1000,), 953),
    (euler50, (), 997651),
    (euler51, (2,), 2),
    (euler51, (3,), 2),
    (euler51, (4,), 2),
    (euler51, (5,), 11),
    (euler51, (6,), 13),
    (euler51, (7,), 56003),
    (euler51, (), 121313),
    (euler52, (2,), 125874),
    (euler52, (), 142857),
    (euler53, (22,), 0),
    (euler53, (23,), 4),
    (euler53, (24,), 11),
    (euler53, (), 4075),
    (euler54, (), 376),
    (euler56, (), 972),
    (euler57, (10,), 1),
    (euler57, (), 153),
    (euler58, (), 26241),
    (euler75, (48, ), 6),
    (euler75, (), 161667),
    (euler62, (3,), 41063625),
    (euler62, (), 127035954683),
    (euler63, (), 49),
    (euler69, (10,), 6),
    (euler69, (), 510510),
    (euler70, (), 8319823),
    (euler72, (8,), 21),
    (euler72, (), 303963552391),
    (euler85, (18,), 6),
    (euler85, (), 2772),
    (euler87, (50,), 4),
    (euler87, (50000000,), 1097343),
    (euler90, (), 1217),
    (euler91_bruteforce, (2,), 14),
    (euler91, (2,), 14),
    (euler91, (), 14234),
    (euler97, (), 8739992577),
    (euler100, (), 756872327473),
    (euler104_, (False, False), 1),
    (euler104_, (False, True), 541),
    (euler104_, (True, False), 2749),
    # (euler104_, (True, True), "?"),
    (euler112_, (90,), 21780),
    # (euler112_, (), "?"),
    (euler113, (6,), 12951),
    (euler113, (10,), 277032),
    (euler113, (), 51161058134250),
    (euler114, (7,), 17),
    (euler114, (), 16475640049),
    (euler115, (3,), 30),
    (euler115, (), 168),
    (euler118, (), 44680),
    (euler119, (2,), 512),
    (euler119, (10,), 614656),
    (euler119, (), 248155780267521),
    (euler120, (), 333082500),
    (euler124, (10, 4), 8),
    (euler124, (10, 6), 9),
    (euler124, (), 21417),
    (euler125, (1000,), 4164),
    (euler125, (), 2906969179),
    (euler127, (1000,), 12523),
    (euler127, (), 18407904),
    (euler131, (100, ), 4),
    (euler141, (100000, ), 124657),
    (euler141, (), 878454337159),
    (euler131, (), 173),
    (euler164, (), 378158756814587),
    (k_max_for_euler183, (11,), 4),
    (k_max_for_euler183, (8,), 3),
    (euler183, (100,), 2438),
    (euler183, (), 48861552),
    (euler191, (4,), 43),
    (euler191, (30,), 1918080160),
    (euler214, (20, 4), 12),
    (euler214, (40000000, 25), 1677366278943),
    (euler215, (9, 3), 8),
    (euler215, (), 806844323190414),
    (euler225, (1,), 27),
    (euler225, (), 2009),
    (euler243, (), 892371480),
    (euler301, (), 2178309),
    (euler387, (4, ), 90619),
    # (euler387, (), 696067597313468),  # Not fast enough - see TODO
    (euler491, (), 194505988824000),
    (euler500_, (4,), 120),
    (euler500_, (5,), 840),
    (euler500_, (6,), 7560),
    (euler500_, (7,), 83160),
    (euler500_, (8,), 1081080),
    # (euler500_, (9,), 17297280),
    # (euler500_, (10,), 294053760),
    # (euler500_, (200,), 482867762),
    # (euler500_, (5000,), 212898634),
    # (euler500_, (10000,), None),
    # (euler500_, (50000,), None),
    # (euler500_, (100000,), None),
    # (euler500_, (), None),
    (euler684, (), 922058210),
    (euler686, (12, 1), 7),
    (euler686, (12, 2), 80),
    (euler686, (123, 45), 12710),
    (euler686, (), 193060223),
    (euler700, (17*1249, 12071), None),
    (euler700, (), None),
    (euler719, (10**4,), 41333),
    (euler719, (10**10,), 499984803177),
    (euler719, (), 128088830547982),
    (euler759, (), None),
    (euler788, (2, ), 18),
    (euler788, (3, ), 270),
    (euler788, (4, ), 603),
    (euler788, (5, ), 8307),
    (euler788, (10, ), 21893256),
    # (euler788, (), 471745499),  # Not fast enough yet
]


def run_tests(tested_func=None):
    """Run tests"""
    nb_fail = 0
    sum_time = 0
    for func, arg, res in tests:
        if tested_func is None or tested_func == func:
            start = timer()
            ret = func(*arg)
            end = timer()
            time = end - start
            sum_time += time
            if ret == res:
                print(func.__name__, "OK", "%.3f" % time)
            else:
                nb_fail += 1
                print(func.__name__, "exp:", res, "got:", ret, "%.3f" % time)
    print("%d failures in %.3f" % (nb_fail, sum_time))
    return nb_fail


def main():
    """Main function"""
    return run_tests()

if __name__ == "__main__":
    exit(main())