|
%I #15 Dec 07 2019 12:18:27
%S 1,5,9,8,4,7,9,9,5,2,16,20,62,8,3,9,28,9,44,95,58,3,2,5,8,8,7,28,0,10,
%T 59,9,7,4,78,6,6,60,6,54,66,9,0,60,9,2,7,0,1,88,0,96,9,0,6,6,0,0,305,
%U 6,4,9,9,94,270,7,9,2,6,93,1,3,5,7,6,9,35,57,9,8,0
%N Whole number sieve of Pi.
%H Manfred Scheucher, <a href="/A247747/b247747.txt">Table of n, a(n) for n = 1..552</a>
%H Manfred Scheucher, <a href="/A247747/a247747.sage.txt">Sage Script</a>
%e Find the first occurrence of 0 (the first whole number) in the digits of Pi (only 35 digits in this illustration):
%e 31415926535897932384626433832795028..., and replace it with a space:
%e 31415926535897932384626433832795 28... Repeat the process with the next whole number, 1:
%e 3 415926535897932384626433832795 28... Then 2:
%e 3 4159 6535897932384626433832795 28... Then 3:
%e 4159 6535897932384626433832795 28... Then 4,5,6,7, etc., until the first occurrence of every counting number is eliminated from the digits of Pi.
%e 1 5 9 8 4 ... Then consolidate gaps between the remaining digits into a single comma:
%e 1,5,9,8,4,7,9,9,5,2,16,20,6,8,3,9, ... to produce the first terms in the whole number sieve of Pi.
%o (Python)
%o def arccot(x, unity):
%o sum = xpower = unity // x
%o n = 3
%o sign = -1
%o while 1:
%o xpower = xpower // (x*x)
%o term = xpower // n
%o if not term:
%o break
%o sum += sign * term
%o sign = -sign
%o n += 2
%o return sum
%o def pi(digits):
%o unity = 10**(digits + 10)
%o pi = 4 * (4*arccot(5, unity) - arccot(239, unity))
%o return pi // 10**10
%o def primes(n):
%o """ Returns a list of primes < n """
%o sieve = [True] * n
%o for i in range(3, int(n**0.5)+1, 2):
%o if sieve[i]:
%o sieve[i*i::2*i]=[False]*((n-i*i-1)/(2*i)+1)
%o return [2] + [i for i in range(3, n, 2) if sieve[i]]
%o a = pi(400)
%o b = range(100000)
%o y = str(a)
%o for x in b:
%o if str(x) in y:
%o y = y.replace(str(x), " ", 1)#replace first occurrence only
%o while " " in y:
%o y = y.replace(" ", " ")#replace long chains of spaces with a single space
%o z = y.split(" ")#split terms into a list
%o z = filter(None, z)#remove null terms
%o f = map(int, z)#convert to integers
%o print(f[0:-1])
%o # Code for A245770 by David Consiglio, Jr., Jan 03 2015
%o # Modified by Manfred Scheucher, Jun 05 2015
%Y Cf. A000796, A245770, A258640, A258481, A257835.
%K nonn,base
%O 1,2
%A _Gil Broussard_, Sep 23 2014
%E Corrected and extended by _Manfred Scheucher_, Jun 05 2015
|
