0249 素数子集和
* * * *
拉格朗日计划
* * * *
素数子集和

记$S=\{2, 3, 5, \ldots, 4999\}$为小于5000的所有素数组成的集合。

在S的所有子集中,有些子集的元素和为素数,求这样的子集的数量的最后十六位数字。

本题难度:



解答

容易验证(比如用Wolfram Alpha)$S$中共有这669个元素,所有元素之和是1548136。

先用筛法生成所有不超过1548136的素数,接下来使用动态规划递推计算。

设$f(i,j)$是在前i个素数组成的集合中和为j的子集的总数,记$p_i$为第i个素数,$t_i$为前i个素数之和,则有 $$f(i,j)=\begin{cases} 1 & j=0, \\ f(i-1,j) & j < p_i, \\ f(i-1,j)+f(i-1,j-p_i) & p_i\le j\le t_{i-1}, \\ f(i-1,j-p_i) & t_{i-1} < j\le t_i. \end{cases}$$ 最后令p取遍1548136以内的素数,计算$f(669,p)$之和并取最后16位即可,结果是$9275262564250418$。

注:计算量较大,代码中打印了进度信息。

target=1548136
d=[0]*target
n=2
while n < target:
    for i in range(n*n,target,n):
        d[i]=d[i]+1
    n=n+1
    while n < target and d[n]>0:
        n=n+1

primeList=[k for k in range(2,target) if d[k]==0]

print "prime list generated with length", len(primeList)

S=primeList[:669]

f=[1,0,1]
t=2

for i in range(1,669):
    t+=S[i]
    f=[f[j] if j < S[i] else (f[j]+f[j-S[i]] if j <=t-S[i] else f[j-S[i]]) for j in range(t+1)]
    if i%7==0:print "progress:",i/7,"percent"

print str(sum(f[j] for j in primeList))[-16:]