题目27：找出为连续数字产生最多质数的二次公式

欧拉计划

Quadratic primes

题目：

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



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

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

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



其中 |n| 表示 n 的绝对值。

例如， |11| = 11， |-4| = 4

对于能够为从 0 开始的连续的 n 产生最多数量的质数的二次公式，找出该公式的系数乘积。

purplenight

愤怒的大头菇

import math
def prime(x):
      if x > 2:
            if x % 2 == 0:
                  return False
            if x == 2:
                  return True
            for i in range(3,int(math.sqrt(x))+1,2):
                  if x % i == 0:
                        return False
            return True
      return False
list2 = []
for i in range(-999,1000):
      for j in range(-999,1000):
            count = 0
            list1 = []
            while True:
                  if prime(count**2+i*count+j):
                        list1.append(count)
                        count += 1
                  else:
                        break
                  
            if list1 not in list2:
                  list2.append(list1)
first = 0
for each in list2:                  #得到list2之后遍历列表找到最长的一项
      if first < len(each):   
            first = len(each)
print(first)                        

得到最长的时候是71个连续的数 从0到70

import math
def prime(x):
      if x > 2:
            if x % 2 == 0:
                  return False
            if x == 2:
                  return True
            for i in range(3,int(math.sqrt(x))+1,2):
                  if x % i == 0:
                        return False
            return True
      return False
list2 = []
for i in range(-999,1000):
      for j in range(-999,1000):
            count = 0
            list1 = []
            while True:
                  if prime(count**2+i*count+j):
                        list1.append(count)
                        if count == 70:   #根据得到的first得到i和j
                              print(i,j)   
                        count += 1
                  else:
                        break
                  
            if list1 not in list2:
                  list2.append(list1)
first = 0
for each in list2:                  #得到list2之后遍历列表找到最长的一项
      if first < len(each):   
            first = len(each)
print(first) 

得到a=-61 b=971
答案是-59231

jerryxjr1220

根据题目的已知条件推断，n最长不可能超过题目中的79，所以只有计算0~79即可。
另外当n=0时，得b必然是要质数。所以列出1000以内的质数供判断。
再来看a，必然是奇数，而且只可能在-b与b之间。
71 -61 971 -59231
[Finished in 2.1s]

def getPrime(n):
        primes = [True]*n
        primes[0],primes[1]=False,False
        for (i, prime) in enumerate(primes):
                if prime:
                        for j in range(i*i,n,i): primes[j] = False
        return [k for (k,trueprime) in enumerate(primes) if trueprime]

def judgePrime(n):
        testlist = getPrime(int(n**0.5+1))
        for each in testlist:
                if n%each == 0:
                        return False
        return True

blist = getPrime(1000)
countmax = 0
for b in blist:
        for a in range(-b,b,2):
                count = 0
                for n in range(80):
                        if n*n+a*n+b <2: break
                        else:
                                if judgePrime(n*n+a*n+b): count += 1
                                else: break
                if countmax < count:
                        countmax, maxa, maxb = count, a, b
print (countmax,maxa,maxb,maxa*maxb)

芒果加黄桃

# encoding:utf-8
# n**2 + a*n + b 找到得到连续素数最多的a、b系数
from time import time
from math import sqrt 
       
def getPrimes(N):
    primes = [True] * N
    primes[0], primes[1] = False, False
    for i, prime in enumerate(primes):
        if prime:
            for j in range(i * i, N, i):
                primes[j] = False
    return [k for k, p in enumerate(primes) if p]
def isPrime(N, prime_list):
    sqr = sqrt(N) + 1
    for t in prime_list:
        if N % t == 0:
            return False 
        if t > sqr:
            break
    return True
def euler027():
    prime_list = getPrimes(1000) 
    count_max = 0
    for b in prime_list:
        for a in range(-b, b, 2):
            count = 0
            for n in range(80):
                t = n ** 2 + a * n + b
                if t < 2:
                    break
                else:
                    if isPrime(t, prime_list):
                        count += 1
                    else:
                        break
            if count_max < count:
                count_max, max_a, max_b = count, a, b
    print(count_max, max_a, max_b)
if __name__ == '__main__':
    start = time()
    euler027()
    print('cost %.6f sec' % (time() - start))

学习了4楼老大的方法  感谢~

渡风

此代码使用matlab编程
Problem27所用时间为16.7542秒
Problem27的答案为-59231

%% Problem27
% 题目27：找出能够带入连续数字得到最多质数的二次公式n²+an+b
%细节：将0带入，得到b必为质数，且大于0
%将1带入，得到a+b+1为质数，a为奇数数，且-b<a<b
function Output=Problem27(Input)
  tic
if nargin==0
Input=1000;
end
b=[];
for ii=Input:-1:3
    if isprime(ii)==1
        Temp=[b,ii];%b为1000以下的质数
        b=Temp;
    end
end
Num=0;
Max=0;
n=0;
Flag=1;
for jj=b
    Temp=jj;
    a=-jj:2:jj;
    for kk=a
        while (Flag==1)
            if n^2+kk*n+jj>0&&isprime(n^2+kk*n+jj)==1
                Num=Num+1;
                n=n+1;
            else
                Flag=0;
            end
        end
        if Num>Max
            Max=Num;
            A=jj;
            B=kk;
        end
        Flag=1;
        Num=0;
        n=0;
    end
end
Output=A*B;
toc
disp('此代码使用matlab编程')
disp(['Problem27所用时间为',num2str(toc),'秒'])
disp(['Problem27的答案为',num2str(Output)])

99592938

def isprime(x):
    if x<=1:return 0
    else:
        for i in range(2,int(x**0.5)+1):
            if not x%i:return 0
        return 1

temp=a=b=0
for i in range(-999,1000):
    for j in range(-999,1000):
        n=0
        f=j
        while isprime(f):
            n+=1
            f=n**2+i*n+j
        if temp<n:
            temp=n
            a=i
            b=j
print(a*b,a,b,temp)

-59231 -61 971 71
>>>

由我们主宰

public class Quadratic{
        public static void main(String[] args){
                int a = -1000,b = -1000,result = 0,n = 0;
                int count = 0,countmax = 0;
                int amax = 0,bmax = 0;
                for(a = -1000;a <= 1000;a ++){
                        for(b = -1000;b <= 1000;b ++){
                                count = 0;
                                for(n = 0;;n ++){
                                        result = n*n+a*n+b;
                                        if(Prime(result)){
                                                count ++;
                                        }
                                        else{
                                                break;
                                        }
                                }
                                if(countmax < count){
                                        amax = a;
                                        bmax = b;
                                        countmax = count-1;
                                }
                        }
                }
                System.out.println(amax+"\t"+bmax+"\t"+countmax+"\t");
                System.out.println("所求最大乘积为"+(amax*bmax));
        }
        public static boolean Prime(int result){
                boolean flag =true;
                for(int i = 2;i < Math.sqrt(Math.abs(result));i ++){
                        if(result % i == 0){
                                flag = false;
                                break;
                        }
                }
                return flag;
        }
}

-61     971     71
所求最大乘积为-59231

王小召

最大连续质数：71
方程：n**2 + -61*n + 971
a*b：59231
用时: 4.5708293 秒

import math
import time


def is_prime(num):
    if num == 2 or num == 3:
        return True
    elif num < 0 or num == 1 or not num % 2 or not num % 3:
        return False
    else:
        for i in range(2, int(math.sqrt(num))+1):
            if not num % i:
                return False
        return True

result = {}
for a in range(-1000, 1000):
    for b in range(-1000, 1000):
        ix = 0
        while is_prime(ix**2 + a*ix + b):
            ix += 1
        result[ix] = (a, b)

print("最多产生连续质数：{} 个\n对应方程：N**2 + {}*N + {}\na*b：{}".format(
    max(result),
    result[max(result)][0], result[max(result)][1],
    abs(result[max(result)][0]*result[max(result)][1])))

print("用时: {} 秒".format(time.process_time()))

永恒的蓝色梦想

from time import process_time as t
t1=t()
from math import sqrt
def isprime(x):
    if x in(2,3):return True
    if x<4 or x%2==0 or(x%6!=1 and x%6!=5):return False
    for i in range(3,int(sqrt(x)+1),2):
        if x%i==0:return False
    return True
l=(0,0,0)
for a in range(-999,1000):
    for b in range(-999,1000):
        for n in range(1,9999):
            if not isprime(n**2+a*n+b):break
        if n>l[2]:l=(a,b,n)
print(f'运行结果：{l[0]*l[1]}\n运行时间：{t()-t1}s')

运行结果：-59231
运行时间：6.640625s

debuggerzh

-59231

Process returned 0 (0x0)   execution time : 0.212 s
Press any key to continue.
参考了4#的思路，然鹅只推出b为奇质数……
只好枚举

#include<iostream>
#include<cmath>
#include<cstdlib>
using namespace std;

const int M = 1000;
const int INF = 1e6;
int cnt = 0;
bool prime[M];
int factor[M];

void euler_sieve(){
  prime[1] = false;
  for (int i = 2;i <= M;i++) prime[i] = true;

  for (int i = 2;i <= M;i++){
    if (prime[i]) factor[++cnt] = i;
    for (int j = 1;j <= cnt && i*factor[j] < M;j++){
      prime[i*factor[j]] = false;
      if (i % factor[j] == 0) break;
    }
  }
}

bool isprime(int x){
  x = abs(x);
  if (x == 0 || x == 1) return false;
  for (int i = 1;factor[i] <= sqrt(x);i++)
    if (x % factor[i] == 0)  return false;
  return true;
}

void tst(){
  for (int i = 0;i < cnt;i++) cout << factor[i] << "*";
  cout << endl;
}

int main(){
  euler_sieve();
  //tst();
  int ab;
  int max_len = 0;

  for (int a = -1000 + 1;a < 1000;a++){
    for (int i = 2;factor[i] < 1000;i++){
      int b = factor[i];

      for (int n = 1;n < 80;n++){
        int fx = n*n + a*n + b;
        if (!isprime(fx)){
          if (n > max_len)  {//cout << a << " " << b << " " << n << endl;
            ab = a*b;  max_len = n;
          }
          break;
        }
      }

      for (int n = 1;n < 80;n++){
        int fx = n*n + a*n - b;
        if (!isprime(fx)){
          if (n > max_len)  {//cout << a << " " << -b << " " << n << endl;
            ab = a*b;  max_len = n;
          }
          break;
        }
      }
    }
  }
  //cout << max_len << " ";
  cout << ab << endl;
}

a1351468657

#include <stdio.h>
#include <math.h>

int check_prime(int n);
int check_prime(int n)
{
        int i, k;
        k = sqrt(n + 1);
        for (i = 2; i <= k; i++)
        {
                if (n % i == 0)
                {
                        return 0;
                }
        }
        return 1;
}

main()
{
        int i, j, n, a, b[200], e = 1, flag;
        int max = 0, x, y;
        b[0] = 2;
        for (i = 3; i < 1000; i++)//根据公式b一定是小于1000的质数
        {
                flag = check_prime(i);
                if (flag == 1)
                {
                        b[e++] = i;
                        
                }
        }
        
        for (a = 0; a < 1000; a++)
        {
                for (i = 0; i < e; i++)
                {
                        n = 1;
                        while (1)
                        {
                                j = n * n - a * n + b[i];
                                if (j < 0)
                                {
                                        break;
                                }
                                flag = check_prime(j);

                                if (flag == 0)
                                {
                                        break;
                                }
                                n++;
                        }
                        if (max < n)
                        {
                                max = n;
                                x = a;
                                y = i;
                        }                
                }
        }
        printf("%d, %d, %d", max, x, b[y]);//x保存最大时a的值，y保存b的当时的下标        
}

ft215378

#找出为连续数字产生最多质数的二次公式
from math import *
from time import *
def is_prime(number):
    if number > 1:
        if number == 2:
            return True
        if number % 2 == 0:
            return False
        for a in range(3,int(sqrt(number) + 1),2):
            if number % a == 0:
                return False
        return True
    return False

start = time()
#当n == 0 时 b一定是质数
b_list = []
for i in range(1000):
    if is_prime(i) and i > 79:
        b_list.append(i)
#print(b_list)
#print(len(b_list))

a = -999
ab_list = []
while True:
    for b in b_list:
        if is_prime(1 + a + b):
            ab_list.append((a, b))   
    a += 1
    if a == 999:
        break
n = 0
index = 0
numlist = [0]
index_list = []
length = len(ab_list)
while index < length:
    if is_prime((n**2) + (n*ab_list[index][0]) + ab_list[index][1]):
        n += 1
    else:
        numlist.append(n)
        index_list.append(index)
        index += 1
        n = 0
        if index == length -1:
            break
index = numlist.index(max(numlist)) - 1 #第n个数是最大的  索引从0开始
#print(ab_list[index])
result = ab_list[index][0] * ab_list[index][1]
end = time()
print("可以连续产生%d个质数" % max(numlist))
print("系数a：%d, 系数b：%d" % (ab_list[index][0], ab_list[index][1]))
print("系数的乘积为：%d" % result)
print("用时：%.2f s" % (end - start))

'''   


#验证

for i in range(71):
    if is_prime((i ** 2) + (i * (-61) + 971)):
        pass
    else:
        print(i)
        print("不对")
        break
    print("bingo")
#'''
      

TLM

python
-59231
Time:0.1665503978729248s
主要应用了b的取值为质数，a为大于-b的奇数
计算是否为奇数，提前计算了较大范围的质数集合，判断是否在其中是非常快的

import time
tt=time.time()
def getzhishu(num):
    l=list(range(1,num+1))
    l[0]=0
    s=set()
    for i in range(num):
        if not l[i]==0:
            s.add(i+1)
            n=2
            while n*(i+1)<=num:
                l[n*(i+1)-1]=0
                n=n+1
    return s
s=getzhishu(79**2+1000*79+1000)
best=(0,0)
for b in getzhishu(1000):
    for a in range(-b+2,1000,2):
        # print((a,b))
        for n in range(80):
            if (n+a)*n+b not in s:
                best=max(best,(n-1,a*b))
                break

print(best[1])
print('Time:{}s'.format(time.time()-tt))

mathtimes

-59231  0.036s
使用 primal 库，Miller-Rabin方法检测质数

extern crate primal;
use primal::*;

fn f(a: i64, b: i64) -> usize {
    for n in 0.. {
        if (a + n) * n + b <= 1 || ((a + n) * n + b > 1 && !is_prime(((a + n) * n + b) as u64)) {
            return n as usize;
        }
    }
    0
}

fn main() {
    let mut max_i = 0;
    let mut max_j = 0;
    let mut max_n = 0;
    for i in -999..999 {
        for j in -999..999 {
            let n = f(i as i64, j as i64);
            if n > max_n {
                max_i = i;
                max_j = j;
                max_n = n;
            }
        }
    }
    println!("b: {}\nc: {}\nn: {}\nans: {}", max_i, max_j, max_n, max_i * max_j);
}