欧拉计划 23

如果一个数的所有真因子之和等于这个数，那么这个数被称为完全数。
例如，28的所有真因子之和为1 + 2 + 4 + 7 + 14 = 28，所以28是一个完全数。

如果一个数的所有真因子之和小于这个数，称其为不足数，如果大于这个数，称其为过剩数。

12是最小的过剩数，1 + 2 + 3 + 4 + 6 = 16。因此最小的能够写成两个过剩数之和的数字是24。

经过分析，可以证明所有大于28123的数字都可以被写成两个过剩数之和。
但是这个上界并不能被进一步缩小，即使我们知道最大的不能表示为两个过剩数之和的数字要比这个上界小。

找出所有不能表示为两个过剩数之和的正整数之和。

def get_divisors(x):
    """ 
    求x所有真因子
    若x = x1 * x2，则x所有真因子 = (x1所有真因子 + x1) 和 (x2所有真因子 + x2) 交叉乘积集合
    """
    x_s = {1}
    x_sqrt = int(pow(x, 0.5))
    for x1 in range(2, x_sqrt + 1):
        if x % x1 == 0:
            x1_s = divisors_dict[x1].copy()
            x1_s.add(x1)
            x2 = x // x1
            x2_s = divisors_dict[x2].copy()
            x2_s.add(x2)
            for i1 in x1_s:
                for i2 in x2_s:
                    x_s.add(i1 * i2)
            x_s.remove(x)
            break
    return x_s


# 求取过剩数
n = 10000
divisors_dict = dict()
for i in range(2, 28123 - 12):
    if i not in divisors_dict:
        divisors_dict[i] = get_divisors(i)
abundant_set = set()
for (k, v) in divisors_dict.items():
    v_sum = sum(v)
    if v_sum > k:
        abundant_set.add(k)
abundant_list = list(abundant_set)
abundant_list.sort()

#sum_2_abundant_set = set()
#for i in range(len(abundant_list)):
#    if abundant_list[i] + abundant_list[i] > 28123:
#        break
#    for j in range(i, len(abundant_list)):
#        if abundant_list[i] + abundant_list[j] > 28123:
#            break
#        sum_2_abundant_set.add(abundant_list[i] + abundant_list[j])
#print(sum(set(range(1, 28124)) - sum_2_abundant_set))
#t1 = time.time()
#print(t1 - t0)

# 当前数 - 过剩数 是否为 过剩数
not_abundant_sum = sum(range(1, 28124))
for i in range(1, 28124):
    for j in abundant_list:
        if j * 2 > i:
            break
        if i - j in abundant_set:
            not_abundant_sum -= i
            break
print(not_abundant_sum)



# 方式2，一次循环，判断正整数是否为两个过剩数之和，及是否为过剩数
abundant_set = set()
abundant_2_not_sum = (1 + 28123) * 28123 / 2
for i in range(2, 28124):
    if len(abundant_set) != 0:
        for a in abundant_set:
            if i - a in abundant_set:
                abundant_2_not_sum -= i
                break
    i_sum = 1
    for j in range(2, int(pow(i, 0.5)) + 1):
        if i % j == 0:
            i_sum += j + i // j
    if int(pow(i, 0.5)) ** 2 == i:
        i_sum -= int(pow(i, 0.5))
    if i_sum > i:
        abundant_set.add(i)
print(abundant_2_not_sum)