login

Year-end appeal: Please make a donation to the OEIS Foundation to support ongoing development and maintenance of the OEIS. We are now in our 61st year, we have over 378,000 sequences, and we’ve reached 11,000 citations (which often say “discovered thanks to the OEIS”).

A245770
Prime sieve of Pi.
7
1, 9, 6, 5, 93, 84626, 3, 2, 502884, 69399, 5105820, 4944, 2, 816406286, 986, 4825342, 70, 9, 480865, 2, 66, 93, 55058, 25, 4081284, 174, 270, 85, 555964462, 4895, 9644288, 7566, 344, 28475648, 8, 65, 201, 45648, 23, 4543, 393, 2602, 412, 724, 660, 55, 74881, 20962, 25, 1715364, 892, 600
OFFSET
1,2
LINKS
Robert G. Wilson v, Table of n, a(n) for n = 1..1029 (first 422 terms from Manfred Scheucher)
Manfred Scheucher, Sage Script (Note: should also run in pure python with a few modifications)
EXAMPLE
Find the first occurrence of prime 2 in the digits of Pi (only 25 digits in this illustration):
3141592653589793238462643..., and replace it with a space:
314159 653589793238462643... Repeat the process with 3:
14159 653589793238462643... Then 5:
141 9 653589793238462643... Then 7:
141 9 653589 93238462643... Then 11, 13, 17, etc., until the first occurrence of every prime is eliminated from the digits of Pi.
1 9 6 5 93 84626 ... Then consolidate gaps between the remaining digits into a single comma:
1,9,6,5,93,84626,... to produce the first terms in the prime sieve of Pi.
PROG
(Python)
def arccot(x, unity):
....sum = xpower = unity // x
....n = 3
....sign = -1
....while 1:
........xpower = xpower // (x*x)
........term = xpower // n
........if not term:
............break
........sum += sign * term
........sign = -sign
........n += 2
....return sum
def pi(digits):
....unity = 10**(digits + 10)
....pi = 4 * (4*arccot(5, unity) - arccot(239, unity))
....return pi // 10**10
def primes(n):
....""" Returns a list of primes < n """
....sieve = [True] * n
....for i in range(3, int(n**0.5)+1, 2):
........if sieve[i]:
............sieve[i*i::2*i]=[False]*((n-i*i-1)/(2*i)+1)
....return [2] + [i for i in range(3, n, 2) if sieve[i]]
a = pi(370)
b = primes(100000000)
y = str(a)
for x in b:
....if str(x) in y:
........y = y.replace(str(x), " ", 1)#replace first occurrence only
........
while " " in y:
....y = y.replace(" ", " ")#replace long chains of spaces with a single space
z = y.split(" ")#split terms into a list
z = filter(None, z)#remove null terms
f = map(int, z)#convert to integers
print(f[0:-1])
# David Consiglio, Jr., Jan 03 2015
CROSSREFS
Cf. A000796, A246022, A248804 (Prime sieve of e), A248831 (Prime sieve of sqrt(2)).
Sequence in context: A144665 A019884 A103814 * A117020 A011459 A354295
KEYWORD
nonn,base
AUTHOR
Gil Broussard, Aug 01 2014
EXTENSIONS
a(14) corrected by Manfred Scheucher, May 25 2015
STATUS
approved