cp's OEIS Frontend

A245770 Prime sieve of Pi.

%I A245770 #38 May 22 2025 10:21:40
%S A245770 1,9,6,5,93,84626,3,2,502884,69399,5105820,4944,2,816406286,986,
%T A245770 4825342,70,9,480865,2,66,93,55058,25,4081284,174,270,85,555964462,
%U A245770 4895,9644288,7566,344,28475648,8,65,201,45648,23,4543,393,2602,412,724,660,55,74881,20962,25,1715364,892,600
%N A245770 Prime sieve of Pi.
%H A245770 Robert G. Wilson v, <a href="/A245770/b245770.txt">Table of n, a(n) for n = 1..1029</a> (first 422 terms from Manfred Scheucher)
%H A245770 Manfred Scheucher, <a href="/A245770/a245770_1.sage.txt">Sage Script</a> (Note: should also run in pure python with a few modifications)
%e A245770 Find the first occurrence of prime 2 in the digits of Pi (only 25 digits in this illustration):
%e A245770 3141592653589793238462643..., and replace it with a space:
%e A245770 314159 653589793238462643...  Repeat the process with 3:
%e A245770 14159 653589793238462643...  Then 5:
%e A245770 141 9 653589793238462643...  Then 7:
%e A245770 141 9 653589 93238462643...  Then 11, 13, 17, etc., until the first occurrence of every prime is eliminated from the digits of Pi.
%e A245770 1   9 6  5   93  84626  ...  Then consolidate gaps between the remaining digits into a single comma:
%e A245770 1,9,6,5,93,84626,...         to produce the first terms in the prime sieve of Pi.
%o A245770 (Python)
%o A245770 def arccot(x, unity):
%o A245770     sum = xpower = unity // x
%o A245770     n = 3
%o A245770     sign = -1
%o A245770     while 1:
%o A245770         xpower = xpower // (x*x)
%o A245770         term = xpower // n
%o A245770         if not term:
%o A245770             break
%o A245770         sum += sign * term
%o A245770         sign = -sign
%o A245770         n += 2
%o A245770     return sum
%o A245770 def pi(digits):
%o A245770     unity = 10**(digits + 10)
%o A245770     pi = 4 * (4*arccot(5, unity) - arccot(239, unity))
%o A245770     return pi // 10**10
%o A245770 def primes(n):
%o A245770     """ Returns  a list of primes < n """
%o A245770     sieve = [True] * n
%o A245770     for i in range(3,int(n**0.5)+1,2):
%o A245770         if sieve[i]:
%o A245770             sieve[i*i::2*i]=[False]*((n-i*i-1)/(2*i)+1)
%o A245770     return [2] + [i for i in range(3,n,2) if sieve[i]]
%o A245770 a = pi(370)
%o A245770 b = primes(100000000)
%o A245770 y = str(a)
%o A245770 for x in b:
%o A245770     if str(x) in y:
%o A245770         y = y.replace(str(x)," ",1)#replace first occurrence only
%o A245770         
%o A245770 while "  " in y:
%o A245770     y = y.replace("  "," ")#replace long chains of spaces with a single space
%o A245770 z = y.split(" ")#split terms into a list
%o A245770 z = filter(None, z)#remove null terms
%o A245770 f = map(int,z)#convert to integers
%o A245770 print(f[0:-1])
%o A245770 # _David Consiglio, Jr._, Jan 03 2015
%Y A245770 Cf. A000796, A246022, A248804 (Prime sieve of e), A248831 (Prime sieve of sqrt(2)).
%K A245770 nonn,base
%O A245770 1,2
%A A245770 _Gil Broussard_, Aug 01 2014
%E A245770 a(14) corrected by _Manfred Scheucher_, May 25 2015