欧拉计划 27

欧拉曾发表过一个著名的二次公式：

n² + n + 41

这个公式对于0到39的连续数字能够产生40个质数。
但是当 n = 40时，402 + 40 + 41 = 40(40 + 1) + 41能够被41整除。
当n = 41时, 41² + 41 + 41显然也能被41整除。

利用计算机，人们发现了一个惊人的公式：n² − 79n + 1601。
这个公式对于n = 0 到 79能够产生80个质数。
这个公式的系数，−79 和1601的乘积是−126479。

考虑如下形式的二次公式：

n² + an + b, 其中|a| < 1000， |b| < 1000

其中|n| 表示 n 的绝对值。
例如， |11| = 11， |−4| = 4
对于能够为从0开始的连续的n产生最多数量的质数的二次公式，找出该公式的系数乘积。

#def get_primes(n):
#    """
#    获取n以下所有质数
#    从2开始，由小到大，依次划去当前数的倍数
#    增加标志list，为0表示划去的数字
#    """
#    flag_list = [1] * n
#    flag_list[0] = 0
#    flag_list[1] = 0
#    n_set = set()
#    for i in range(n):
#        if flag_list[i]:
#            n_set.add(i)
#            for j in range(i * 2, n, i):
#                flag_list[j] = 0
#    return n_set


def get_primes(n):
    """
    利用set差集，获取n以下所有质数
    从2开始，由小到大，依次划去当前数的倍数    
    """
    n_set = set(range(3, n, 2))
    n_set.add(2)
    for i in range(3, n, 2):
        if i in n_set:
            i_set = set(range(i * 2, n, i))
            n_set -= i_set
    return n_set    


def is_prime(x):
    """ 判断x是否是质数 """
    if x < 2:
        return False
    if x == 2:
        return True
    x_sqrt = int(pow(x, 0.5))
    l = [2]
    l.extend(range(3, x_sqrt + 1, 2))
    for i in l:
        if x % i == 0:
            return False
    return True


fn = lambda n, a, b: n ** 2 + a * n + b
ab_max = 1000
n_max = 1
a_prod_b = 0

# n = 0 => b为1000内质数
# n = 1 => a > 2 - b
for b in get_primes(ab_max):
    for a in range(2 - b, ab_max):
        n = 1
#        # 长度不够，继续
#        if not is_prime(fn(n_max, a, b)):
#            continue
        while 1:
            if not is_prime(fn(n, a, b)):
                break
            n += 1
        if n > n_max:
            n_max = n
            a_prod_b = a * b
print(a_prod_b)