N:= 2*10^7: # for terms <= N
P:= select(isprime, [2, seq(i, i=3..N/2, 2)]):
nP:= nops(P):
SP:= 0:
for i from 1 while P[i]^2 <= N do
m:= ListTools:BinaryPlace(P, N/P[i]);
SP:= SP, op(P[i]*P[i..m]);
od:
SP:= sort([SP]):
SS:= ListTools:PartialSums(SP):
V:= Vector(6):
SI:= Vector(6):
II:= Vector(6, 1):
for i from 1 to 6 do SI[i]:= SS[i+1]SS[1] od:
count:= 1: V[1]:= 4: m:= 4: im:= {1}:
while count < 6 do
for j in im do
II[j]:= II[j]+1;
SI[j]:= SS[II[j]+j]  SS[II[j]];
od;
m:= min(SI);
im:= select(j > SI[j] = m, {$1..20});
for k from 1 to 20 while {$1..k} subset im do
if V[k] = 0 then V[k]:= m; count:= count+1 fi
od;
od:
convert(V, list);
