%I #39 Jul 25 2017 17:12:46
%S 2,3,5,7,13,19,23,31,43,47,73,83,109,113,181,199,283,293,313,317,463,
%T 467,503,509,523,619,661,683,691,887,1063,1069,1103,1109,1123,1129,
%U 1303,1307,1321,1327,1613,1621,1627,1637,1669,1789
%N Popular primes: primes prime(n) such that prime(n-k)+prime(n+k) >= 2*prime(n) for all k = 1,2,...n-2.
%C McNew says that a prime p is "popular" on an interval [2, k] if no prime occurs more frequently than p as the greatest prime factor (gpf, A006530) of the integers in that interval. - _N. J. A. Sloane_, Jul 25 2017
%C The first two primes, 2 and 3, are tested against an empty set of k, and we include them, defining such a test to have a positive outcome.
%H Chai Wah Wu, <a href="/A124661/b124661.txt">Table of n, a(n) for n = 1..10000</a>
%H Nathan McNew, <a href="http://arxiv.org/abs/1504.05985">Popular values of the largest prime divisor function</a>, arXiv:1504.05985 [math.NT], 2015.
%H Nathan McNew, <a href="http://dx.doi.org/10.1080/10586458.2016.115518">The Most Frequent Values of the Largest Prime Divisor Function</a>, Exper. Math., 2017, Vol. 26, No. 2, 210-224.
%H C. Pomerance, <a href="http://dx.doi.org/10.1090/S0025-5718-1979-0514836-7">The prime number graph</a>, Math. Comp. 33 (1979) 399--408. - _Nathan McNew_, Apr 04 2014
%e prime(11)=31 is in the sequence because prime(10)+prime(12) = 66, prime(9)+prime(13) = 64,..., prime(2)+prime(20) = 74 are all >= 62 = 2*31.
%e prime(10) = 29 is not in the sequence because prime(9)+prime(11) = 54 for example is smaller than 58 = 2*29.
%t Select[Prime@ Range@ 300, Function[{p, n}, NoneTrue[Range[n - 2], Prime[n - #] + Prime[n + #] < 2 p &]] @@ {#, PrimePi@ #} &] (* _Michael De Vlieger_, Jul 25 2017 *)
%o (PARI) isok(p) = {n = primepi(p); for (k=1, n-2, if (prime(n-k) + prime(n+k) < 2*p, return (0));); return (1);}
%o lista(nn) = {for(n=1, nn, if (isok(prime(n)), print1(prime(n), ", ");););} \\ _Michel Marcus_, Nov 03 2013
%o (Python)
%o def a124661(end):
%o ....a124661_list=[2, 3, 5, 7]
%o ....primes=[2,3]
%o ....add=True
%o ....for num in range(4, end*3):
%o ........prime=False
%o ........length=len(primes)
%o ........for y in range(0,length):
%o ............if num % primes[y]!=0:
%o ................prime=True
%o ............else:
%o ................prime=False
%o ................break
%o ........if (prime):
%o ............primes.append(num)
%o ....for x in range(4, int(len(primes)/2)):
%o ........for k in range(1, x-2):
%o ............if (primes[x-k]+primes[x+k]>=primes[x]*2):
%o ................add=True
%o ............else:
%o ................add=False
%o ................break
%o ........if (add):
%o ............if (primes[x]>end):
%o ................break
%o ............else:
%o ................a124661_list.append(primes[x])
%o ....return a124661_list
%o # _Conner L. Delahanty_, Apr 19 2014
%o (Python)
%o from sympy import prime
%o A124661_list = []
%o for n in range(1,10**6):
%o p = prime(n)
%o for k in range(1,n-1):
%o if prime(n-k)+prime(n+k) < 2*p:
%o break
%o else:
%o A124661_list.append(p) # _Chai Wah Wu_, Jul 25 2017
%Y Cf. A006530, A051635.
%K nonn,easy
%O 1,1
%A _Artur Jasinski_, Dec 23 2006
%E Sequence extended by _R. J. Mathar_, Mar 28 2010
|