Machin-like formula

From HandWiki
Short description: Formulas for pi


In 1706, William Jones published a contribution from John Machin of the first 100 decimal digits of the circle constant π = 3.14159..., at that time a record for accuracy, along with the unexplained and surprising formula used to calculate them,[1][2]

π=(1654239)13(165342393)+15(165542395).

This formula is an expanded form of the equation now called Machin's formula,

π4=4arctan15arctan1239,

where the arctangents of 15 and 1239 have been expanded using the arctangent series,

arctanx=xx33+x55x77+.

When applied directly to find 14π=arctan1, the arctangent series converges extremely slowly, requiring five billion terms to obtain 10 correct decimal digits. Machin's formula is dramatically more practical, needing only six terms to obtain 10 correct digits.

Several other well-known mathematicians immediately set to work making sense of Machin's formula, developing their own variants and extensions, now in general called Machin-like formulas. These have the form

c0π4=n=1Ncnarctananbn

 

 

 

 

(1)

where c0 is a positive integer, cn are signed non-zero integers, and an and bn are positive integers such that an<bn.

Derivation

The angle addition formula for arctangent asserts that

arctana1b1+arctana2b2=arctana1b2+a2b1b1b2a1a2,

 

 

 

 

(3)

if π2<arctana1b1+arctana2b2<π2. All of the Machin-like formulas can be derived by repeated application of equation 3. As an example, we show the derivation of Machin's original formula. One has: 2arctan15=arctan15+arctan15=arctan15+155511=arctan1024=arctan512, and consequently 4arctan15=2arctan15+2arctan15=arctan512+arctan512=arctan512+512121255=arctan120119. Therefore also 4arctan15π4=4arctan15arctan11=4arctan15+arctan11=arctan120119+arctan11=arctan1201+(1)1191191120(1)=arctan1239, and so finally π4=4arctan15arctan1239.

An insightful way to visualize equation 3 is to picture what happens when two complex numbers are multiplied together:

(b1+a1i)(b2+a2i)
=b1b2+a2b1i+a1b2ia1a2

=(b1b2a1a2)+(a1b2+a2b1)i

 

 

 

 

(4)

The angle associated with a complex number (bn+ani) is given by:

arctananbn

Thus, in equation 4, the angle associated with the product is:

arctana1b2+a2b1b1b2a1a2

Note that this is the same expression as occurs in equation 3. Thus equation 3 can be interpreted as saying that multiplying two complex numbers means adding their associated angles (see multiplication of complex numbers).

The expression:

cnarctananbn

is the angle associated with:

(bn+ani)cn

Equation 1 can be re-written as:

k(1+i)c0=n=1N(bn+ani)cn

Here k is an arbitrary constant that accounts for the difference in magnitude between the vectors on the two sides of the equation. The magnitudes can be ignored, only the angles are significant.

Using complex numbers

Other formulas may be generated using complex numbers.[3] For example, the angle of a complex number (a+bi) is given by arctanba and, when one multiplies complex numbers, one adds their angles. If a=b then arctanba is 45 degrees or π4 radians. This means that if the real part and complex part are equal then the arctangent will equal π4. Since the arctangent of one has a very slow convergence rate if we find two complex numbers that when multiplied will result in the same real and imaginary part we will have a Machin-like formula. An example is (2+i) and (3+i). If we multiply these out we will get (5+5i). Therefore, arctan12+arctan13=π4.

In order to use complex numbers to show that π4=4arctan15arctan1239, it is first necessary to understand that raising a complex number to a real power k implies multiplying its anomaly (angle) by k, and that the anomaly of the product of two complex numbers is equal to the sum of their anomalies. Since it can by shown, by doing the calculation, that (5+i)4(239i)=(1+i)22134, i.e. that the real and imaginary parts of both sides are equal, and since that equality is equivalent to: 4arctan15arctan1239=π4, the latter equality is also demonstrated.

Lehmer's measure

One of the most important parameters that characterize computational efficiency of a Machin-like formula is the Lehmer's measure, defined as[4][5]

λ=n=1N1log10(bn/an).

In order to obtain the Lehmer's measure as small as possible, it is necessary to decrease the ratio of positive integers an/bn in the arctangent arguments and to minimize the number of the terms in the Machin-like formula. Nowadays at an=1 the smallest known Lehmer's measure is λ1.51244 due to H. Chien-Lih (1997),[6] whose Machin-like formula is shown below. It is very common in the Machin-like formulas when all numerators an=1.

Two-term formulas

In the special case where the numerator an=1, there are exactly four solutions having only two terms.[7][8] All four were found by John Machin in 1705–1706, but only one of them became widely known when it was published in William Jones's book Synopsis Palmariorum Matheseos, so the other three are often attributed to other mathematicians. These are

Euler's 1737 (known to Machin 1706):[9][10]

π4=arctan12+arctan13

Hermann's 1706 (known to Machin 1706):[11][10]

π4=2arctan12arctan17

Hutton's or Vega's (known to Machin 1706):[8][10]

π4=2arctan13+arctan17

and Machin's 1706:[1][10]

π4=4arctan15arctan1239 .

In the general case, where the value of a numerator an is not restricted, there are infinitely many other solutions. For example:

π4=22arctan128+arctan174450748218032836685456512798646395734210062276153190241239

or

π4=22arctan24478873121+17arctan68560169049993

 

 

 

 

(5)

Example

The adjacent diagram demonstrates the relationship between the arctangents and their areas. From the diagram, we have the following:

area(PON)=area(MOF)=π×MOF2π=MEF=arctan12area(POM)=area(NOF)=arctan13area(POF)=π4=area(PON)+area(NOF)=arctan12+arctan13area(MON)=arctan17area(PON)=arctan12=area(POM)+area(MON)=arctan13+arctan17,

a relation which can also be found by means of
the following calculation within the complex numbers

(3+i)(7+i)=211+(3+7)i=10(2+i).

More terms

The 2002 record for digits of π, 1,241,100,000,000, was obtained by Yasumasa Kanada of Tokyo University. The calculation was performed on a 64-node Hitachi supercomputer with 1 terabyte of main memory, performing 2 trillion operations per second. The following two equations were both used:

π4=12arctan149+32arctan1575arctan1239+12arctan1110443
Kikuo Takano (1982).
π4=44arctan157+7arctan123912arctan1682+24arctan112943
F. C. M. Størmer (1896).

Two equations are used so that one can check they both give the same result; it is helpful if the equations used to cross-check the result reuse some of the arctangent arguments (note the reuse of 57 and 239 above), so that the process can be simplified by only computing them once, but not all of them, in order to preserve their independence.

Machin-like formulas for π can be constructed by finding a set of m integers bn,n=1..m, where all the prime factorisations of bn2+1, taken together, use a number of distinct primes m, and then using either linear algebra or the LLL basis-reduction algorithm to construct linear combinations of arctangents of 1bn. For example, in the Størmer formula above, we have

572+1=25313
2392+1=2134
6822+1=53612
129432+1=25413361

so four expressions whose factors are powers of only the four primes 2, 5, 13 and 61.

In 1993 Jörg Uwe Arndt[12] found the 11-term formula:

π4=36462arctan1390112+135908arctan1485298+274509arctan168398239581arctan11984933+178477arctan12478328114569arctan13449051146571arctan118975991+61914arctan12270927469044arctan12420814489431arctan120122958243938arctan12189376182

using the set of 11 primes {2,5,13,17,29,37,53,61,89,97,101}.

Another formula where 10 of the arctan-arguments are the same as above has been discovered by Hwang Chien-Lih (黃見利) (2004), so it is easier to check they both give the same result:

π4=36462arctan151387+26522arctan1485298+19275arctan16839823119arctan119849333833arctan124783285183arctan1344905137185arctan11897599111010arctan122709274+3880arctan12420814416507arctan12012295827476arctan12189376182

It can be observed that these formulas reuse all the same arctangents after the first one. They are constructed by looking for numbers where b2+1 is divisible only by primes less than 102.

The most efficient currently known Machin-like formula for computing π is:

π4=183arctan1239+32arctan1102368arctan15832+12arctan111044312arctan14841182100arctan16826318
(Hwang Chien-Lih, 1997)

where the set of primes is {2,5,13,229,457,1201}.

A further refinement is to use "Todd's Process", as described in;[5] this leads to results such as

π4=183arctan1239+32arctan1102368arctan15832+12arctan1113021100arctan1682631812arctan133366019650+12arctan143599522992503626068
(Hwang Chien-Lih, 2003)

where the large prime 834312889110521 divides the bn2+1 of the last two indices.
M. Wetherfield found 2004

π4=83arctan1107+17arctan1171022arctan110369724arctan1251348944arctan118280007883+12arctan17939642926390344818+22arctan13054211727257704725384731479018.

On Pi Day 2024, Matt Parker along with 400 volunteers used the following formula to hand calculate π:

π4=1587arctan12852+295arctan14193+593arctan14246+359arctan139307+481arctan155603+625arctan1211050708arctan1390112

It was the biggest hand calculation of π in a century. [13]

More methods

There are further methods to derive Machin-like formulas for π with reciprocals of integers. One is given by the following formula:[14]

π4=2k1arctan1Ak+m=1Marctan1Bk,m+arctan1Bk,M+1,

where

a0:=0

and recursively

ak:=2+ak1,Ak:=ak2ak1

and

Bk,1:=2(Ak+iAki)2k1ii

and recursively

Bk,m:=1+Bk,m1Bk,m1Bk,m1Bk,m1.

E.g., for k=4 and M=5 we get:

π4=8arctan110arctan184arctan121342arctan1991268848arctan1193018008592515208050arctan1197967899896401851763240424238758988350338arctan1117573868168175352930277752844194126767991915008537018836932014293678271636885792397

This is verified by the following MuPAD code:

z:=(10+I)^8*(84-I)*(21342-I)*(991268848-I)*(193018008592515208050-I)\
  *(197967899896401851763240424238758988350338-I)\
  *(117573868168175352930277752844194126767991915008537018836932014293678271636885792397-I):
Re(z)-Im(z)
0

meaning

z:=(10+i)8(84i)(21342i)(991268848i)(193018008592515208050i)(197967899896401851763240424238758988350338i)(117573868168175352930277752844194126767991915008537018836932014293678271636885792397i)=(1+i)(z).

Efficiency

For large computations of π, the binary splitting algorithm can be used to compute the arctangents much, much more quickly than by adding the terms in the Taylor series naively one at a time. In practical implementations such as y-cruncher, there is a relatively large constant overhead per term plus a time proportional to 1/logbn, and a point of diminishing returns appears beyond three or four arctangent terms in the sum; this is why the supercomputer calculation above used only a four-term version.

It is not the goal of this section to estimate the actual run time of any given algorithm. Instead, the intention is merely to devise a relative metric by which two algorithms can be compared against each other.

Let Nd be the number of digits to which π is to be calculated.

Let Nt be the number of terms in the Taylor series (see equation 2).

Let un be the amount of time spent on each digit (for each term in the Taylor series).

The Taylor series will converge when:

((bnan)2)Nt=10Nd

Thus:

Nt=Ndln102lnbnan

For the first term in the Taylor series, all Nd digits must be processed. In the last term of the Taylor series, however, there's only one digit remaining to be processed. In all of the intervening terms, the number of digits to be processed can be approximated by linear interpolation. Thus the total is given by:

NdNt2

The run time is given by:

time=unNdNt2

Combining equations, the run time is given by:

time=unNd2ln104lnbnan=kunlnbnan

Where k is a constant that combines all of the other constants. Since this is a relative metric, the value of k can be ignored.

The total time, across all the terms of equation 1, is given by:

time=n=1Nunlnbnan

un cannot be modelled accurately without detailed knowledge of the specific software. Regardless, we present one possible model.

The software spends most of its time evaluating the Taylor series from equation 2. The primary loop can be summarized in the following pseudo code:

1: term *= an2
2: term /= bn2
3: tmp = term / (2*n+1)
4: sum += tmp

In this particular model, it is assumed that each of these steps takes approximately the same amount of time. Depending on the software used, this may be a very good approximation or it may be a poor one.

The unit of time is defined such that one step of the pseudo code corresponds to one unit. To execute the loop, in its entirety, requires four units of time. un is defined to be four.

Note, however, that if an is equal to one, then step one can be skipped. The loop only takes three units of time. un is defined to be three.

As an example, consider the equation:

π4=44arctan7468414967113+139arctan123912arctan2013815351991

 

 

 

 

(6)

The following table shows the estimated time for each of the terms:

an bn bnan lnbnan un time
74684 14967113 200.41 5.3003 4 0.75467
1 239 239.00 5.4765 3 0.54780
20138 15351991 762.34 6.6364 4 0.60274

The total time is 0.75467 + 0.54780 + 0.60274 = 1.9052

Compare this with equation 5. The following table shows the estimated time for each of the terms:

an bn bnan lnbnan un time
24478 873121 35.670 3.5743 4 1.1191
685601 69049993 100.71 4.6123 4 0.8672

The total time is 1.1191 + 0.8672 = 1.9863

The conclusion, based on this particular model, is that equation 6 is slightly faster than equation 5, regardless of the fact that equation 6 has more terms. This result is typical of the general trend. The dominant factor is the ratio between an and bn. In order to achieve a high ratio, it is necessary to add additional terms. Often, there is a net savings in time.

References

  1. 1.0 1.1 Jones, William (1706). Synopsis Palmariorum Matheseos. London: J. Wale. pp. 243, 263. https://archive.org/details/SynopsisPalmariorumMatheseosOrANewIntroductionToTheMathematics/page/n283/. "There are various other ways of finding the Lengths, or Areas of particular Curve Lines or Planes, which may very much facilitate the Practice; as for instance, in the Circle, the Diameter is to Circumference as 1 to
    165423913165342393+15165542395,&c.=
    3.14159, &c. = π. This Series (among others for the same purpose, and drawn from the same Principle) I receiv'd from the Excellent Analyst, and my much Esteem'd Friend Mr. John Machin; and by means thereof, Van Ceulen's Number, or that in Art. 64.38. may be Examin'd with all desireable Ease and Dispatch."
     

    Reprinted in Smith, David Eugene (1929). "William Jones: The First Use of π for the Circle Ratio". A Source Book in Mathematics. McGraw–Hill. pp. 346–347. https://archive.org/details/sourcebookinmath1929smit/page/346/. 

  2. Beckmann, Petr (1971). A History Of Pi. USA: The Golem Press. p. 102. ISBN 0-88029-418-3. https://archive.org/details/historyofpisymbo00beck/page/102. 
  3. Størmer, Carl (1897). "Sur l'application de la théorie des nombres entiers complexes a la solution en nombres rationnels x1 x2xnc1 c2cnk de l'équation: c1arc tgx1+c2arc tgx2++cnarc tgxn=kπ4". Archiv for Mathematik og Naturvidenskab 19 (3): 1–95. https://archive.org/details/archivformathema1918961897oslo/page/n122/. 
  4. Lehmer, Derrick Henry (1938). "On Arccotangent Relations for π". American Mathematical Monthly 45 (10): 657–664. doi:10.2307/2302434. 
  5. 5.0 5.1 Wetherfield, Michael (2016). "The Enhancement of Machin's Formula by Todd's Process". The Mathematical Gazette 80 (488): 333–344. doi:10.2307/3619567. 
  6. Chien-Lih, Hwang. "More Machin-Type Identities". The Mathematical Gazette 81 (490). 
  7. Størmer, Carl (1896). "Solution complète en nombres entiers m, n, x, y et k de l'équation marc tg1x+narc tg1y=kπ4.". Skrifter Udgivne Af Videnskabsselskabet I Christiania 1895 (11): 1–21. https://archive.org/details/skrifterudgivnea1895chri/page/n750/. 
  8. 8.0 8.1 Størmer, Carl (1899). "Solution complète en nombres entiers de l'équation marctang1x+narctang1y=kπ4" (in fr). Bulletin de la Société Mathématique de France 27: 160–170. doi:10.24033/bsmf.603. https://archive.org/details/bulletinsocit27soci/page/n179/. 
  9. Euler, Leonhard (1744). "De variis modis circuli quadraturam numeris proxime exprimendi". Commentarii Academiae Scientiarum Petropolitanae 9: 222–236. E 74. https://archive.org/details/commentariiacade09impe/page/222/. 
  10. 10.0 10.1 10.2 10.3 Tweddle, Ian (1991). "John Machin and Robert Simson on Inverse-tangent Series for π". Archive for History of Exact Sciences 42 (1): 1–14. doi:10.1007/BF00384331. 
  11. Letter from Jakob Hermann to Gottfried Leibniz, 21 August 1706. Published in Gerhardt, C.I., ed (1859). "XXII. Hermann an Leibniz.". Leibnizens mathematische Schriften. 4. H.W. Schmidt. pp. 302–304. https://archive.org/details/leibnizensmathe06leibgoog/page/302/. 
  12. Jörg Uwe Arndt: "Matters Computational" section 32.5.2, page 637.
  13. The biggest hand calculation in a century! [Pi Day 2024]. Retrieved 2024-04-02 – via www.youtube.com.
  14. https://arxiv.org/pdf/2108.07718.pdf (2021)