www.pi314.net


History
Mathematicians
All formulas
Num. approx.
Softwares
Misc. math.
Digits
Poetry
Papers/videos
Delirium !
 Pi-Day
Images
Music
Links
Bibliography



Boris Gourévitch
The world of Pi - V2.57
modif. 13/04/2013



Pi-Day in
Pi Day Countdown
Google
Home Version history Guestbook Who I am Some pictures (fr) Acknowledgements Last modifications Contact

Cette page en français This page in English


Stanley Rabinowitz

The spigot algorithm
A. Sale - D. Saada - S. Rabinowitz



The Principle


So, what is this algorithm?
Well, let's remember that we know pi in its decimal form, which is in base 10:

Pi=3.141592653589793238462643383279...
this can also be written as :

This is known as Horner's form; it allows us, among other things, to reduce the number of multiplications when evaluating a polynomial. We can see here that it is base 10, and that the base's step is constant, meaning that we have 10 at each decimal place (digit if base different from 10). We note this base [1/10,1/10,1/10...]
Now, let's consider the series discovered by Euler (See page about him for the demonstration) and let's write it in Horner's form:


Now this is interesting, if we compare it with the previous expression, we can see that we are considering a base with a changing step [1,1/3,2/5,3/7,4/9...]. And in this base, Pi is [2;2,2,2,2...]. So in this base, Pi is one of the simpliest numbers that exists !
We know Pi 's digits in this base, so to compute Pi 's decimal places in base 10 one by one, one just needs to build an algorithm that changes it to base 10, which is precisely the principle of the spigot algorithm.

A historical overview of this method

I might as well tell you that I haven't found much information about the above mathematicians! We can't always find what we want on the web, and if someone doesn't have a personal page, it limits what one can find out about him.
Initially, it was A. Sale who had the idea of this method in 1968 and used it to compute e.
Then D. Saada suggested applying it to Pi in 1988 and so did S. Rabinowitz in 1991. The latter is quite well known, he's a confirmed hacker, and even a mathematician who has published quite a few articles (He is a Mac friend by the way). He studied some of the finesses of this algorithm with his collegue Stan Wagon in 1995 in an article in the Mathematical American Monthly.

Around :

First of all, a little computation concerning how much memory it will take. In Horner's form, we can see that the n/(2n+1) step of the base is slightly under 1/2 every time. The exact value 1/2 would bring us to consider base 2. For the conversion from base 2 to base 10, we will need roughly Log2(10n)=3,32n digits. We can consider it is about the same value for our base with a changing step to base 10. So if we want to get n décimal places, we will have to consider 3.32n boxes. To get 4 décimal places (including the 3 before the decimal point), we build a memory table of about 14 columns. In fact, 13 will be enough here. The translation of this algorithm would bring us to consider a table like this:


A

r

1

2

3

4

5

6

7

8

9

10

11

12

B

=

 

3

5

7

9

11

13

15

17

19

21

23

25

Initialization   2 2 2 2 2 2 2 2 2 2 2 2 2
                             
*10

20

20

20

20

20

20

20

20

20

20

20

20

20

carried over

10

12

12

12

10

12

7

8

9

0

0

0

0

sum

3

30

32

32

32

30

32

27

28

29

20

20

20

20

remainder

0

2

2

4

3

10

1

13

12

1

20

20

20

                             
*10

0

20

20

40

30

100

10

130

120

10

200

200

200

carried over

13

20

33

40

65

48

98

88

72

150

132

96

0

sum

1

13

40

53

80

95

148

108

218

192

160

332

296

200

remainder

3

1

3

3

5

5

4

8

5

8

17

20

0

                             
*10

30

10

30

30

50

50

40

80

50

80

170

200

0

carried over

11

24

30

40

40

42

63

64

90

120

88

0

0

sum

4

41

34

60

70

90

92

103

144

140

200

258

200

0

remainder

1

1

0

0

0

4

12

9

4

10

6

16

0

                             
*10

10

10

0

0

0

40

120

90

40

100

60

160

0

carried over

4

2

9

24

55

84

63

48

72

60

66

0

0

sum

1

14

12

9

24

55

124

183

138

112

160

126

160

0

remainder

4

0

4

3

1

3

1

3

10

8

0

22

0



And now let's explain how and why it works !
We recognise in the two first lines the numerator and the denominator of the steps of the changing step base. In the third line, we find Pi 's expression in that base. And we fill the last column of the remainder line with 0 's.
The algorithm to convert from right to left in the table is as follows :
Formally, the general principle consists in multiplying the digit by 10 and in computing the remainder in an equivalent of a Euclidian division of this number by the step of the changing step base. A number carried over will then appear at each calculation and comes from the same phenomenon which appears when one multiplies 53 by 8; firstly one multiplies 3 by 8, then one carries over 2 and adds it to the multiplication of 8 by 5 which is the multiplication of the next digit.


Filling in the red column for example :

1) First we fill in the *10 line by multiplying the previous line by 10. No problem so far !
2) We place in sum the sum of the *10 line with the carried over line, which is:

20+9=29
3) We do the Euclidian division of the sum by the number in line B of the same column :

29=17*1+12
4) We place the 12 remainder in the remainder line. Then we multiply the 1 quotient with line A and put the result in the next column's carried over line :

1*8=8
The 9 remainder of the red column comes from the previous column's computation. And so here we pass on the 8 remainder to the next sum's computation.

By repeating the process with the other columns, we fill in the first table and we get 30 as the last sum. But as we multiplied everything by 10, we take 3 as the first digit of Pi. That's it!

One little remark however : We can consider that a number in the last column is right if it is not followed by a 9.
When the remainder in the r column is higher than 100, we can obtain (but it is very rare...) 10 after the last decimal place that we take for Pi. We then take this digit plus 1 as digit for Pi, It's as simple as that!

Other formulas ?

Well I think any rational series should be all right as long as in it's Horner form, we only find integrers as expression of Pi in the changing step base (for the previous series, we only had 2's). Particularly, Ramanujan rational series should be able to be used and give a better algorithm.
Gosper's page's formula is also valid :


As the fraction of two terms of the series is always smaller than 1/13, we will need Log13(10n)=0.9 boxes to get a digit. So if we consider 6 boxes, we can hope to get 6*1/0.9=6.6 decimal places, so 6 or even 7 decimal places at best. We can see that this is respected perfectly in the table because we even obtain 7 décimal places :


A

r

1

6

15

28

45

B

=

 

60

168

330

546

816

Initialization   3 8 13 18 23 28
               
*10

30

80

130

180

230

280

carried over

1

0

0

0

0

0

sum

3

31

80

130

180

230

280

remainder

1

20

130

180

230

280

               
*10

10

200

1300

1800

2300

2800

carried over

4

48

75

112

135

0

sum

1

14

248

1375

1912

2435

2800

remainder

4

8

31

262

251

352

               
*10

40

80

310

2620

2510

3520

carried over

41

12

120

112

180

0

sum

4

41

92

430

2732

2690

3520

remainder

1

32

94

92

506

256

               
*10

10

320

940

920

5060

2560

carried over

5

30

45

252

135

0

sum

1

15

350

985

1172

5195

2560

remainder

5

50

145

182

281

112

               
*10

50

500

1450

1820

2810

1120

carried over

9

54

75

140

45

0

sum

5

59

554

1525

1960

2855

1120

remainder

9

14

13

310

125

304

               
*10

90

140

130

3100

1250

3040

carried over

2

6

135

56

135

0

sum

9

92

146

265

3156

1385

3040

remainder

2

26

97

186

293

592

               
*10

20

260

970

1860

2930

5920

carried over

4

36

90

140

315

0

sum

2

24

296

1060

2000

3245

5920

remainder

4

56

52

20

515

208


back to home page