security/nss/lib/freebl/mpi/doc/pi.txt@6474c204b198 (annotated)
security/nss/lib/freebl/mpi/doc/pi.txt
Wed, 31 Dec 2014 06:09:35 +0100
- author
- Michael Schloh von Bennewitz <michael@schloh.com>
- date
- Wed, 31 Dec 2014 06:09:35 +0100
- changeset 0
- 6474c204b198
- permissions
- -rw-r--r--
Cloned upstream origin tor-browser at tor-browser-31.3.0esr-4.5-1-build1
revision ID fc1c9ff7c1b2defdbc039f12214767608f46423f for hacking purpose.
| michael@0 | 1 | This file describes how pi is computed by the program in 'pi.c' (see |
| michael@0 | 2 | the utils subdirectory). |
| michael@0 | 3 | |
| michael@0 | 4 | Basically, we use Machin's formula, which is what everyone in the |
| michael@0 | 5 | world uses as a simple method for computing approximations to pi. |
| michael@0 | 6 | This works for up to a few thousand digits without too much effort. |
| michael@0 | 7 | Beyond that, though, it gets too slow. |
| michael@0 | 8 | |
| michael@0 | 9 | Machin's formula states: |
| michael@0 | 10 | |
| michael@0 | 11 | pi := 16 * arctan(1/5) - 4 * arctan(1/239) |
| michael@0 | 12 | |
| michael@0 | 13 | We compute this in integer arithmetic by first multiplying everything |
| michael@0 | 14 | through by 10^d, where 'd' is the number of digits of pi we wanted to |
| michael@0 | 15 | compute. It turns out, the last few digits will be wrong, but the |
| michael@0 | 16 | number that are wrong is usually very small (ordinarly only 2-3). |
| michael@0 | 17 | Having done this, we compute the arctan() function using the formula: |
| michael@0 | 18 | |
| michael@0 | 19 | 1 1 1 1 1 |
| michael@0 | 20 | arctan(1/x) := --- - ----- + ----- - ----- + ----- - ... |
| michael@0 | 21 | x 3 x^3 5 x^5 7 x^7 9 x^9 |
| michael@0 | 22 | |
| michael@0 | 23 | This is done iteratively by computing the first term manually, and |
| michael@0 | 24 | then iteratively dividing x^2 and k, where k = 3, 5, 7, ... out of the |
| michael@0 | 25 | current figure. This is then added to (or subtracted from) a running |
| michael@0 | 26 | sum, as appropriate. The iteration continues until we overflow our |
| michael@0 | 27 | available precision and the current figure goes to zero under integer |
| michael@0 | 28 | division. At that point, we're finished. |
| michael@0 | 29 | |
| michael@0 | 30 | Actually, we get a couple extra bits of precision out of the fact that |
| michael@0 | 31 | we know we're computing y * arctan(1/x), by setting up the multiplier |
| michael@0 | 32 | as: |
| michael@0 | 33 | |
| michael@0 | 34 | y * 10^d |
| michael@0 | 35 | |
| michael@0 | 36 | ... instead of just 10^d. There is also a bit of cleverness in how |
| michael@0 | 37 | the loop is constructed, to avoid special-casing the first term. |
| michael@0 | 38 | Check out the code for arctan() in 'pi.c', if you are interested in |
| michael@0 | 39 | seeing how it is set up. |
| michael@0 | 40 | |
| michael@0 | 41 | Thanks to Jason P. for this algorithm, which I assembled from notes |
| michael@0 | 42 | and programs found on his cool "Pile of Pi Programs" page, at: |
| michael@0 | 43 | |
| michael@0 | 44 | http://www.isr.umd.edu/~jasonp/pipage.html |
| michael@0 | 45 | |
| michael@0 | 46 | Thanks also to Henrik Johansson <Henrik.Johansson@Nexus.Comm.SE>, from |
| michael@0 | 47 | whose pi program I borrowed the clever idea of pre-multiplying by x in |
| michael@0 | 48 | order to avoid a special case on the loop iteration. |
| michael@0 | 49 | |
| michael@0 | 50 | ------------------------------------------------------------------ |
| michael@0 | 51 | This Source Code Form is subject to the terms of the Mozilla Public |
| michael@0 | 52 | # License, v. 2.0. If a copy of the MPL was not distributed with this |
| michael@0 | 53 | # file, You can obtain one at http://mozilla.org/MPL/2.0/. |
