题目60：找出一个5个质数的集合，该集合中每两个质数相连接都能产生一个素数
Prime pair sets
The primes 3, 7, 109, and 673, are quite remarkable. By taking any two primes and concatenating them in any order the result will always be prime. For example, taking 7 and 109, both 7109 and 1097 are prime. The sum of these four primes, 792, represents the lowest sum for a set of four primes with this property.
Find the lowest sum for a set of five primes for which any two primes concatenate to produce another prime.
题目：
质数 3, 7, 109, 和 673 是值得注意的。将其中任意两个质数以任何顺序相连接产生的结果都是质数。例如，取 7 和 109，连接而成的 7109 和 1097 都是质数。这四个质数的和是 792，这也是满足这个性质的四个质数集合的最小总和。
找出满足这个性质的五个质数的集合中，集合数之和最小的。算出这个最小的和。
Mark 本帖最后由 jerryxjr1220 于 2016-10-14 09:34 编辑
Result = 26033, in 36.276000 s
#!/usr/bin/env python
"""
q060
Algorithm (by iwek7, from forum):
1. Generate prime list up to N. Remove 2 and 5 from list since they are trivial.
Call the list prime_list.
2. Create an empty list all_list, to contain lists of prime pairs.
3. For each prime, ci, in prime_list, do
a. create an empty new_list to contain newly found lists of prime pairs
b. for each list l in all_list, do
i. Test if ALL primes in l form prime pair with ci.
ii. If it does,
1. Create a copy of list l, call it cl
2. Append ci to this cl
3. If length of cl >= LEN, solution found.
4. Append cl to new_list
c. Append all lists in new_list to all_list
d. Append to all_list
Optimisation:
1. Create a is_prime_table to store the flag if a given prime pairs are prime or not. This
is to avoid primality test for a given prime pair.
2. Use a faster prime generation algorithm.
3. Use a faster primality test algorithm, rather than generating list of primes and test
if a prime pair is in the list or not.
"""
import time
import numpy as np
def getPrime(n):
primes = *n
primes, primes = False, False
for (i,prime) in enumerate(primes):
if prime:
for j in range(i*i,n,i):
primes = False
return primes
def is_prime(n):
factor = 3
if n == 2:
return True
if (n%2 == 0):
return False
while factor * factor <= n:
if n % factor == 0:
return False
factor += 2
return True
def test_cat(prime_set, ci, i, lut):
for j, cj in prime_set:
if lut == 0:
cat1 = int(str(ci) + str(cj))
cat2 = int(str(cj) + str(ci))
isprime = is_prime(cat1) and is_prime(cat2)
lut = 1 if isprime else -1
lut = 1 if isprime else -1
if lut == -1:
return False
return True
def solve_q060(prime, LEN):
all_list = list()
for i, ci in enumerate(prime):
new_list = list()
for l in all_list:
if test_cat(l, ci, i, is_prime_table):
new_cand = list(l)
new_cand.append((i, ci))
if len(new_cand) >= LEN:
return new_cand
new_list.append(new_cand)
all_list.extend(new_list)
all_list.append([(i, ci)])
return None
primes = getPrime(100000)
UPPER = 10000
LEN = 5
tic = time.time()
prime_cand = []
pl = getPrime(UPPER)
for (m,trueprime) in enumerate(pl):
if trueprime:
prime_cand.append(m)
prime_cand.remove(2)
prime_cand.remove(5)
prime_len = len(prime_cand)
res_list = []
is_prime_table = np.zeros((prime_len, prime_len))
res_list = solve_q060(prime_cand, LEN)
res_list = for c in res_list]
toc = time.time()
if res_list:
res = sum(res_list)
print "Result = %d, in %f s" % (res, toc - tic)
else:
print "Not found, try increase prime numbers." 递归算法，20秒不到可以出结果：
13
5197
5701
6733
8389
def is_prime(n):
factor = 2
while factor * factor <= n:
if n % factor == 0:
return False
factor += 1
return True
def is_valid_pair(n, m):
return is_prime(int(str(n) + str(m))) and is_prime(int(str(m) + str(n)))
def check_comb(nb_chosen, indices):
if nb_chosen == NB_TO_CHOOSE:
for i in indices:
print(prime_list)
exit()
start = -1
if nb_chosen > 0:
start = indices
for i in range(start + 1, nb_primes):
valid = True
for j in range(nb_chosen):
if not is_valid_pair(prime_list, prime_list]):
valid = False
break
if valid:
indices = i
check_comb(nb_chosen + 1, indices)
UPPER_BOUND = 10 ** 4
NB_TO_CHOOSE = 5
prime_list = []
prime = * UPPER_BOUND
for i in range(2, UPPER_BOUND):
if prime:
prime_list.append(i)
for j in range(i + i, UPPER_BOUND, i):
prime = False
nb_primes = len(prime_list)
indices = [-1] * NB_TO_CHOOSE
check_comb(0, indices) 另外还有一种字典法，不过用时也不短，关键生成字典的时间比较长。
方法多多 # encoding:utf-8
# 寻找特殊质数3,7,109,673已任何形式组合比如377371091097 都是质数
# 寻找满足条件的和最小的五个质数，满足该条件
from time import time
from math import sqrt
def getPrimes(N):
primes = * N
primes, primes = False, False
for i, prime in enumerate(primes):
if prime:
for j in range(i * i, N, i):
primes = False
return
def isPrime(n, l_pr):
# 不考虑1
if n < 4:
return True
if n % 2 == 0 or n % 3 == 0 :
return False
t = int(sqrt(n))
for i in l_pr:
if i > t:
return True
if n % i == 0:
return False
return True
def findNums(l_pr, l_primes, N, l_result):
if N == 1 and len(l_result) < 1:
l_result.append(min(l_primes))
else:
if len(l_result) < N:
findNums(l_pr, l_primes, N - 1, l_result)
for i in :
flag = True
for j in l_result:
if not (isPrime(int(str(j) + str(i)), l_pr)
and isPrime(int(str(i) + str(j)), l_pr)):
flag = False
break
if flag:
l_result.append(i)
return
def euler060(N=10000):
l_primes = getPrimes(N)
l_result =
l_pr =
findNums(l_pr, l_primes, 5, l_result)
while len(l_result) < 5:
l_primes =
l_result.remove(max(l_result))
findNums(l_pr, l_primes, 5, l_result)
print(l_result)
if __name__ == '__main__':
start = time()
euler060()
print('cost %.6f sec' % (time() - start))
cost 11.987199 sec 用的matlab
结果是
13 5197 5701 6733 8389
时间已过 3.885805 秒。
>>
第二组的话
7 1237 2341 12409 18433
13 5197 5701 6733 8389
时间已过 16.375198 秒。
>> 是谁说 mathematica代码好理解的 下面这段代码，能理解算我输答案26033
pp := FromDigits
Total[FindClique[
Graph[#[] <-> #[] & /@
Select[Prime@Range@PrimePi@10000~Subsets~{2},
PrimeQ[#[]~pp~#[]] && PrimeQ[#[]~pp~#[]] &]]][[
1]]] 感觉很麻烦，但还是出来的不算太慢def prime(n):
if n==2 or n==3:
return True
else:
if n%2==0:
return False
else:
for i in range(3,int(n**0.5)+1,2):
if n%i == 0:
return False
return True
def lj(n,m):#连接两个数字
return prime(int(str(n)+str(m))) and prime(int(str(m)+str(n)))
l1 = list(i for i in range(2,10000) if prime(i))
count = 0
for i in l1:
for j in l1:
if lj(i,j):
for k in l1:
if lj(i,k) and lj(j,k):
for l in l1:
if lj(i,l) and lj(j,l) and lj(k,l):
for m in l1:
if lj(i,m) and lj(j,m) and lj(k,m) and lj(l,m):
count += 1
print(i,j,k,l,m)
if count==1:
exit() 本帖最后由 guosl 于 2022-1-9 22:23 编辑
/*
答案：26033
耗时：30.2571秒（单核）
5.08559秒（六核）
*/
#include <iostream>
#include <vector>
#include <cstring>
#include <cstdlib>
#include <cmath>
#include <omp.h>
using namespace std;
const int N = 30000;//搜索范围一定要大于那个最小和才行，否则从理论上来讲不能得到全局最小和的
const int M = 1000000;
const int INF = 0x7fffffff;
char cp;
vector<int> vp;
bool chkbp(long long n)//判别数n是否是一个素数
{
if ((n & 1) == 0)
return false;
int d = (int)sqrt((double)n);
for (int i = 1; i < (int)vp.size() && vp <= d; ++i)
{
if (n%vp == 0)
return false;
}
return true;
}
bool chk(int p1, int p2)//检查p1p2与p2p1是否是素数
{
char str1, str2;
_itoa_s(p1, str1, 10);
_itoa_s(p2, str2, 10);
char str;
strcpy_s(str, str1);
strcat_s(str, str2);
long long n = atoll(str);
if (n < M)
{
if (cp == 0)
return false;
}
else
{
if (!chkbp(n))
return false;
}
strcpy_s(str, str2);
strcat_s(str, str1);
n = atoll(str);
if (n < M)
{
if (cp == 0)
return false;
}
return chkbp(n);
}
int main(void)
{
double t = omp_get_wtime();
memset(cp + 2, 1, (M + 2) * sizeof(char));
for (int i = 2; i <= (int)sqrt((double)M); ++i)
{
if (cp == 0)
continue;
for (int j = i * i; j <= M; j += i)
cp = 0;
}
vp.push_back(2);
for (int i = 3; i <= N; i += 2)
{
if (cp == 1)
vp.push_back(i);
}
volatile int nMinSum = INF, k = 1;
volatile bool bContinue = true;
#pragma omp parallel shared(k,bContinue,vp,nMinSum)
while (bContinue)
{
if (k >= nMinSum || k >= (int)vp.size())
{
bContinue = false;
#pragma omp flush(bContinue)
break;
}
int i0;
#pragma omp critical(c1)
{
i0 = k;
++k;
}
for (int i1 = i0 + 1; i1 < (int)vp.size(); ++i1)
{
if (vp + vp >= nMinSum)
break;
if (!chk(vp, vp))
continue;
for (int i2 = i1 + 1; bContinue && i2 < (int)vp.size(); ++i2)
{
if (vp + vp + vp >= nMinSum)
break;
if (!(chk(vp, vp) && chk(vp, vp)))
continue;
for (int i3 = i2 + 1; bContinue && i3 < (int)vp.size(); ++i3)
{
if (vp + vp + vp + vp >= nMinSum)
break;
if (!(chk(vp, vp) && chk(vp, vp) && chk(vp, vp)))
continue;
for (int i4 = i3 + 1; bContinue && i4 < (int)vp.size(); ++i4)
{
if (vp + vp + vp + vp + vp >= nMinSum)
break;
if (!(chk(vp, vp) && chk(vp, vp) && chk(vp, vp) && chk(vp, vp)))
continue;
int nS = vp + vp + vp + vp + vp;
if (nMinSum > nS)
#pragma omp critical(c2)
nMinSum = nS;
}
}
}
}
}
t = omp_get_wtime() - t;
cout << nMinSum << endl << t << endl;
return 0;
} import time as t
import math
start = t.perf_counter()
def check_prime(a_num):
is_prime = True
for prime_i in range(2, int(math.sqrt(a_num) + 1)):
if a_num % prime_i == 0:
is_prime = False
return is_prime
def is_exchangeable(prime_list):
length_list = len(prime_list)
for each_prime1 in range(length_list - 1):
for each_prime2 in range(each_prime1 + 1, length_list):
new_num_1 = int(str(prime_list) + str(prime_list))
new_num_2 = int(str(prime_list) + str(prime_list))
if not check_prime(new_num_1) or not check_prime(new_num_2):
return False
return True
def p60():
primes = []
prime_range = 10000
for numbers in range(2, prime_range):
if check_prime(numbers):
primes.append(numbers)
primes.remove(2)
primes.remove(5)
test_range = len(primes)
for prime1 in range(test_range):
for prime2 in range(prime1 + 1, test_range):
cur_plist = , primes]
if is_exchangeable(cur_plist):
for prime3 in range(prime2 + 1, test_range):
cur_plist = , primes, primes]
if is_exchangeable(cur_plist):
for prime4 in range(prime3 + 1, test_range):
cur_plist = , primes, primes, primes]
if is_exchangeable(cur_plist):
for prime5 in range(prime4 + 1, test_range):
cur_plist = , primes, primes,
primes, primes]
if is_exchangeable(cur_plist):
return cur_plist
answer = p60()
print(answer)
print(sum(answer))
print("It costs %f s" % (t.perf_counter() - start))
26033
It costs 51.725813 s
页:
[1]