题目58：考察螺旋网格中对角线上质数的数量

欧拉计划

Spiral primes

Starting with 1 and spiralling anticlockwise in the following way, a square spiral with side length 7 is formed.



It is interesting to note that the odd squares lie along the bottom right diagonal, but what is more interesting is that 8 out of the 13 numbers lying along both diagonals are prime; that is, a ratio of 8/13 ≈ 62%.

If one complete new layer is wrapped around the spiral above, a square spiral with side length 9 will be formed. If this process is continued, what is the side length of the square spiral for which the ratio of primes along both diagonals first falls below 10%?

题目：

从 1 开始逆时针旋转，可以得到一个边长为 7 的螺旋正方形。



有趣的是奇数的平方数都处于右下对角线上。更有趣的是，对角线上的 13 个数字中有 8 个质数，其百分比是 8/13 ≈ 62%。

如果在上面的螺旋正方形上再加一层，可以得到一个边长为 9 的螺旋正方形。如果这个过程继续，到螺旋正方形的边长为多少时，对角线上的质数百分比第一次降到 10% 以下？

迷雾少年

边长为3的正方形的对角线质数有3个 质数百分比为0.333333
边长为5的正方形的对角线质数有5个 质数百分比为0.2
边长为7的正方形的对角线质数有8个 质数百分比为0.163265
边长为9的正方形的对角线质数有9个 质数百分比为0.111111
边长为11的正方形的对角线质数有10个 质数百分比为0.0826446
边长为13的正方形的对角线质数有11个 质数百分比为0.0650888
边长为15的正方形的对角线质数有13个 质数百分比为0.0577778
边长为17的正方形的对角线质数有15个 质数百分比为0.0519031
边长为19的正方形的对角线质数有16个 质数百分比为0.0443213
边长为21的正方形的对角线质数有18个 质数百分比为0.0408163
边长为23的正方形的对角线质数有19个 质数百分比为0.0359168
边长为25的正方形的对角线质数有21个 质数百分比为0.0336
边长为27的正方形的对角线质数有22个 质数百分比为0.0301783
边长为29的正方形的对角线质数有23个 质数百分比为0.0273484
边长为31的正方形的对角线质数有23个 质数百分比为0.0239334
边长为33的正方形的对角线质数有23个 质数百分比为0.0211203
边长为35的正方形的对角线质数有24个 质数百分比为0.0195918
边长为37的正方形的对角线质数有25个 质数百分比为0.0182615
边长为39的正方形的对角线质数有26个 质数百分比为0.017094
边长为41的正方形的对角线质数有27个 质数百分比为0.0160619
边长为43的正方形的对角线质数有28个 质数百分比为0.0151433
边长为45的正方形的对角线质数有28个 质数百分比为0.0138272
边长为47的正方形的对角线质数有28个 质数百分比为0.0126754
边长为49的正方形的对角线质数有28个 质数百分比为0.0116618
边长为51的正方形的对角线质数有29个 质数百分比为0.0111496
边长为53的正方形的对角线质数有29个 质数百分比为0.010324
边长为55的正方形的对角线质数有31个 质数百分比为0.0102479
边长为57的正方形的对角线质数有32个 质数百分比为0.00984918
边长为59的正方形的对角线质数有33个 质数百分比为0.00948003
边长为61的正方形的对角线质数有34个 质数百分比为0.00913733
边长为63的正方形的对角线质数有35个 质数百分比为0.00881834
边长为65的正方形的对角线质数有35个 质数百分比为0.00828402
边长为67的正方形的对角线质数有37个 质数百分比为0.00824237
边长为69的正方形的对角线质数有37个 质数百分比为0.00777148
边长为71的正方形的对角线质数有38个 质数百分比为0.00753819
边长为73的正方形的对角线质数有39个 质数百分比为0.00731845
边长为75的正方形的对角线质数有40个 质数百分比为0.00711111
边长为77的正方形的对角线质数有41个 质数百分比为0.00691516
边长为79的正方形的对角线质数有43个 质数百分比为0.00688992
边长为81的正方形的对角线质数有44个 质数百分比为0.00670629
边长为83的正方形的对角线质数有44个 质数百分比为0.00638699
边长为85的正方形的对角线质数有45个 质数百分比为0.00622837
边长为87的正方形的对角线质数有45个 质数百分比为0.0059453
边长为89的正方形的对角线质数有45个 质数百分比为0.0056811
边长为91的正方形的对角线质数有48个 质数百分比为0.0057964
边长为93的正方形的对角线质数有48个 质数百分比为0.00554977
边长为95的正方形的对角线质数有49个 质数百分比为0.00542936
边长为97的正方形的对角线质数有49个 质数百分比为0.00520778
边长为99的正方形的对角线质数有49个 质数百分比为0.00499949


辣鸡代码还调试了很久，，，，

#include "myLib/myLib.h"
#pragma comment(lib,"myLib//myLib.lib")
using namespace std;


//加入是n行n列 正方形  那么右下角的数字就是n*n   就是质数1-n*n


//建立一个n*n的矩阵 返回其中对角线质数个数
//建立的过程是个递归  从外到里  
//假设最大的变长1000
#define MAXC 1000
int arr[MAXC][MAXC]={0};

int hang = 7;
//建立
//m是第几列了
//startnum是左上角开始的数字

void jianli(int n,int m,int startnum)
{
        if(1==startnum) 
        {
                arr[n][n] = 1;
                return ;
        }
        int temp=startnum;
        for (int z=0;z<m;z++)
        {
                        arr[n][n+z] = temp--;
        }
        temp=startnum;
        for (int z=0;z<m;z++)
        {
                arr[n+z][n] = temp++;
        }
        temp=startnum+m-1;
        for (int z=0;z<m;z++)
        {
                arr[n+m-1][n+z] = temp++;
        }
        temp=arr[n][n+m-1];
        for (int z=0;z<m-1;z++)
        {
                arr[n+z][n+m-1] = temp--;
        }
        jianli(n+1,m-2,startnum-4*m+8);
}

//建立表后找出有多少个对角线质数

int gethow()
{
        int sum=0;
        for (int i = 0;i<hang;i++)
        {
                
                //std::cout<<arr[i][i]<<" "<<endl;

                if(iszhishu(arr[i][i]))
                        sum++;
        }
                int lie  = hang;
                for (int i = 0;i<hang;i++)
                {
                        //std::cout<<arr[i][hang-i-1]<<" "<<endl;

                        if(iszhishu(arr[i][hang-i-1]))
                                sum++;
                }
        return sum;
}
/*
void clear()
{
        hang = 0;
        for (int n=0;n<MAXC;n++)
        {
                for (int nm=0;nm<MAXC;nm++)
                {
                        arr[n][nm] = 0;
                }
        }
}*/
int main()
{
        //3行的正方形 开始的数字是5
        //3 5  3*3 - 4
        //5 17  5*5 -8
        //7 37  7*7 -12
        //开始数字= 边长平方 - 2*(变长-1)
        //表都建立了 接下来就easy了
         
        //std::cout<<gethow();
        for (int i = 3;i<=100;i+=2)
        {
                double baifen=0;
                hang = i;
                jianli(0,i,i*i-2*(i-1));
                std::cout<<"边长为"<<i<<"的正方形"<<"的对角线质数有"<<gethow()<<"个 ";
                baifen = double(gethow() ) /  ((double)i*i);
                std::cout<<"质数百分比为"<<baifen<<endl;
                //clear();
        }
        return 0;
}

jerryxjr1220

迷雾少年 发表于 2016-8-31 12:45
边长为3的正方形的对角线质数有3个 质数百分比为0.333333
边长为5的正方形的对角线质数有5个 质数百分比为 ...

题目理解错了，这个百分比不是对角线质数占所有数的百分比，是占对角线个数的百分比
正确的应该是26241

def isPrime(n):
        if n == 1: return False
        if n == 2 or n == 3: return True
        if n%2 == 0: return False
        for i in range(3,int(n**0.5+1),2):
                if n%i == 0:
                        return False
        return True

crosslist = [0]*100000
for n in range(3,100000,2):
        cross = 2*n-1
        crossprime = 0 
        if isPrime(n*n): crossprime += 1
        if isPrime(n*n-n+1): crossprime += 1
        if isPrime(n*n-2*n+2): crossprime += 1
        if isPrime(n*n-3*n+3): crossprime += 1
        crosslist[n] = crosslist[n-2]+crossprime
        if (float(crosslist[n]) / cross) < 0.1:
                print (n)
                exit()

愤怒的大头菇

import math
def Prime(x):
    if x >1:
        if x==2:
            return True
        if x%2==0:
            return False
        for i in range(3,int(math.sqrt(x)+1),2):
            if x%i==0:
                return False
        return True
    return False
list1=[1];total = 1;count = 0;floor = 1;i = 1
while True:
    i += 2
    floor += 1
    if Prime(i**2):
        count += 1
    if Prime(i**2-(i-1)):
        count += 1
    if Prime(i**2-(i-1)*2):
        count += 1
    if Prime(i**2-(i-1)*3):
        count += 1
    total += 4
    if count /total <0.1:
        break
print(i

结果：26241

芒果加黄桃

# encoding:utf-8
# 计算对角线质数比例
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, 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 euler058(N=100000):
    n_chprimes, step, num = 0, 0, 1
    l_primes = getPrimes(N)
    while True:
        step += 2
        for i in range(1, 5):
            num += step
            flag = isPrime(num, l_primes)
            if flag:
                n_chprimes += 1
        if n_chprimes / (2 * step + 1) < 0.10:
            print(step + 1, n_chprimes / (2 * step + 1))
            break
if __name__ == '__main__':
    start = time() 
    euler058()
    print('cost %.6f sec' % (time() - start))

26241 0.09999809454850327
cost 2.233223 sec

飘飞的白杨

def isprime(n):
    '''判断是否是素数'''
    if n <= 1: return False
    if n == 2 or n == 3: return True
    if n % 2 == 0: return False
    for i in range(3, int(n**0.5)+1, 2):
        if n % i == 0: return False
    return True
num = 0
for i in range(3, 100000, 2):
    if isprime(i*i-(i-1)): num += 1
    if isprime(i*i-2*(i-1)): num += 1
    if isprime(i*i-3*(i-1)): num += 1
    if num/(i*2-1) < 0.1: break
print(i)

改

def primes(N):
    '''素数表'''
    l =  int(N**0.5) + 1
    prime = [False, False] + [True] * N
    for i in range(2, l):
        if prime[i]:
            for j in range(2, N//i + 1):
                if prime[j*i]:
                    prime[j*i] = False
    return tuple(i for i in range(len(prime)) if prime[i])
def isprime(n):
    '''判断是否是素数'''
    if n <= 1: return False
    if n == 2: return True
    for i in primes_tuple:
        if i*i > n: return True
        if n % i == 0: return False
    for i in range(primes_tuple[-1]+2, int(n**0.5)+1, 2):
        if n % i == 0: return False
    return True
X, num = 100000, 0
primes_tuple = primes(X)
for i in range(3, X, 2):
    if isprime(i*i-(i-1)): num += 1
    if isprime(i*i-2*(i-1)): num += 1
    if isprime(i*i-3*(i-1)): num += 1
    if num/(i*2-1) < 0.1: break
print(i)

najin

用的matlab
结果是
       26241

时间已过 4.395559 秒。
>>

刚开始写了函数直接生成螺旋矩阵，算不出来；
后来换了解法，但是判断素数又耗时太长，最后用了这个办法

hk1057

import time


def time_it():
        start=time.time()
        fun()
        end=time.time()
        print('cost %.6f Secs'%(end-start))


def is_prime(num):
        if num <= 1:
                return False
        if num == 2:
                return True
        if num % 2 == 0:
                return False
        for i in range(3,int(num**0.5)+1,2):
                if num % i == 0:
                        return False
        else:
                return True
       
       
def fun():
        start=3
        list1=[1]
        list2=[]
        while True:
                for x in range(4):
                        list1.append(start**2-(start-1)*x)
                       
                for each in list1[-4:]:
                        if is_prime(each):
                                list2.append(each)
               
                if len(list2)/len(list1) < 0.1:
                        print('第一个质数百分比<10%%的螺旋正方形的边长为%d'%start)
                        break
                else:
                        start+=2
                       

if __name__=='__main__':
        time_it()
'''
第一个质数百分比<10%的螺旋正方形的边长为26241
cost 2.991065 Secs
'''

王小召

结果：26241
用时：6.7392432 秒

import time

def is_prime(num):
    if not num % 2 or not num % 3:
        return False
    else:
        for i in range(2, int(num**0.5) + 1, 1):
            if not num % i:
                return False
        return True


def cal():
    fenzi = 8
    l = 7
    while fenzi/(2*l - 1) > 0.1:
        l += 2
        count = 0
        if is_prime(l ** 2):
            count += 1
        if is_prime(l ** 2 - (l - 1)):
            count += 1
        if is_prime(l ** 2 - 2 * (l - 1)):
            count += 1
        if is_prime(l ** 2 - 3 * (l - 1)):
            count += 1
        fenzi += count
    else:
        return l
print("结果：{}\n用时：{} 秒".format(cal(), time.process_time()))

debuggerzh

26241

Process returned 0 (0x0)   execution time : 0.470 s
Press any key to continue.
欧拉筛，找规律，线性递推

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

const int M = 1e5;
int ct = 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[++ct] = i;
    for (int j = 1;j <= ct && i*factor[j] < M;j++){
      prime[i*factor[j]] = false;
      if (i % factor[j] == 0) break;
    }
  }
}

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

int main(){
  euler_sieve();
  int prime_num = 0;

  for (int n = 3; ;n+=2){
    int tot = 2*n - 1;
    int ini = n*n - 3*n + 3;

    for (int i = 0;i < 3;i++){
      bool b = ini < M ? prime[ini] : isprime(ini);
      if (b) prime_num++;

      ini += n-1;
    }//免去对平方数的判断
    if ( (double)prime_num / tot < 0.1 )  {cout << n << endl;  break;}
  }

  return 0;
}

4444567

26241 0.09999809454850327

from time import time
def prime(n):
    if n<= 1:
        return False
    else:
        if n==2 or n==3:
            return True
        if n>3:
            for i in range(3,int(n**0.5)+1,2):
                if n%i == 0:
                    return False
            return True
t = time()
b = 1     #边长,每次+2
djx = 1   #对角线上数字个数，每次+4
c = 1     #初始数字1
d = 0    #递增数
l1= []   #处于对角线的质数，四次为一个周期，2，4，6递增
while True:
    a = 1  #周期指针
    b += 2
    d += 2
    djx += 4
    while a<5:
        c += d
        if prime(c):
            l1.append(c)
        a += 1
    if len(l1)/djx < 0.1:
        print(b, len(l1)/djx)
        print('\n',time()-t)
        break
    

6.351067543029785