The Madhava series is slow and boring, says … Madhava
यत्सङ्ख्ययात्र हरणे कृते निवृत्ता हृतिस्तु जामितया
तस्या ऊर्ध्वगता या समसङ्ख्या तद्दलं गुणोऽन्ते स्यात् ||
तद्वर्गो रूपयुतो हारो व्यासाब्धिघातात् प्राग्वत्
ताभ्यामाप्तं स्वमृणे कृते धने क्षेप एव करणीयः ||
लभ्धः परिधिः सूक्ष्मो बहुकृत्वो हरणतोऽतिसूक्ष्मश्च ||
When you stop division out of boredom, remember the last divisor
take the next even number, and halve it. That is your numerator.
The square of the numerator, plus one is your divisor. Multiply by four times the diameter as before.
The correction you thus calculate, add to the result if the last operation you did was subtraction, and vice versa.
The corrected diameter you get is accurate, much more so that by a lot of division (ie: by taking several more terms!)
What exactly does this mean?
We are back with व्यासे वारिधिनिहते, the Madhava Pi (or more accurately circumference) series. We saw in an earlier article how it was demonstrated to be correct by Yuktibhāṣā. To jog the reader's memory, the verse is :
व्यासे वारिधिनिहते रूपहृते व्याससागराभिहते ।
त्रिशरादिविषमसंख्याभक्तमृणं स्वं पृथक्क्रमात् कुर्यात् ॥
which is equivalent to
\( C = 4.d – \frac{4.d}{3} + \frac{4.d}{5} – \frac{4.d}{7} + … \)
\( \frac{\pi}{4} = 1 – \frac{1}{3} + \frac{1}{5} – \frac{1}{7} \)
The one problem with this series is that it converges extremely slowly. This is precisely why the verse says निवृत्ता हृतिस्तु जामितया (“when you stop dividing out of boredom”). Remember that all this calculation is happening by hand. How many terms would we calculate? 10? 20? 50?
3.14159265358979
These are the first 15 digits of Pi – this is our target.
Madhava's series converges really really slowly towards this number. Setting the diameter to 1, we add more and more terms to the series to calculate the circumference as
| Terms | Result |
| 10 | 3.0418396189294024 |
| 20 | 3.0916238066678385 |
| 30 | 3.1082685666989462 |
| 40 | 3.1165965567938323 |
| 50 | 3.1215946525910105 |
| 100 | 3.131592903558553 |
50 terms is a lot to calculate by hand, and we haven't got anywhere close to the answer we want. Even after a 100 terms, Aryabhata's 3.1416 is still much better. This series as is, is a theoretical curiosity, and is not much use for practical calculation. Or is it?
(We have conveniently ignored that Madhava himself would have calculated in fractions, not decimals, and would've probably set d = 100000000000, but that doesn't change our point)
| Terms | Result | Error | 1/Terms |
| 10 | 3.0418396189294024 | 0.09975303466039076 | 0.1 |
| 20 | 3.0916238066678385 | 0.04996884692195458 | 0.05 |
| 30 | 3.1082685666989462 | 0.03332408689084687 | 0.03333.. |
| 40 | 3.1165965567938323 | 0.02499609679596082 | 0.025 |
| 50 | 3.1215946525910105 | 0.019998000998782572 | 0.02 |
| 100 | 3.131592903558553 | 0.009999750031240318 | 0.01 |
We notice that the difference between Pi and our approximation is almost the same, but a little less than (and grows increasingly close to) the reciprocal of the number of terms we have calculated. If only we could add a correction term …
संस्कारः – The correction term.
Let's start with the presumption that we want a correction term, to be added to the end of the series, depending on the number of terms in the series. The correction term must
- Be dependent on the number of terms already calculated (not a constant!)
- Approximate the rest of the series
- Converge faster than the series itself.
Let's assume we have computed n terms in the series. The last denominator we computed was \(p = 2n-1\). We are looking for a correction term of the form \(1/a_p\) to add to the series, such that
\( C \approx 4.d – \frac{4.d}{3} + \frac{4.d}{5} – \frac{4.d}{7} + … (-1)^{\frac{p-1}{2}}.\frac{4.d}{p} + (-1)^{\frac{p+1}{2}}.\frac{4.d}{a_p} \)
As n (and p) grow larger, we want the series to converge faster with the correction term added. ie, we want:
\( 4.d – \frac{4.d}{3} + \frac{4.d}{5} – \frac{4.d}{7} + … (-1)^{\frac{p-3}{2}}.\frac{4.d}{p-2} + (-1)^{\frac{p-1}{2}}.\frac{4.d}{a_{p-1}} \approx 4.d – \frac{4.d}{3} + \frac{4.d}{5} – \frac{4.d}{7} + … (-1)^{\frac{p-1}{2}}.\frac{4.d}{p} + (-1)^{\frac{p+1}{2}}.\frac{4.d}{a_p} \) which, after accounting for signs, simplifies to:
\( \frac{1}{a_{p-1}} + \frac{1}{a_{p}} \approx \frac{1}{p}\)
We define the स्थौल्य (error term) \( E_p = \frac{1}{a_{p-1}} + \frac{1}{a_{p}} – \frac{1}{p} \), which we want to keep as small as we can, and want it to go to zero with n much faster than \(\frac{1}{n}\) does.
First order correction term
Our first attempt at a correction term is \( a_p = 2(p+1) = 4n \). This is easy to calculate. Whatever your last denominator was, add one to it, and double your result to get the denominator of your correction term. This is a good start, as we can see :
| Terms | Result with first order correction |
| 10 | 3.141839618929402 |
| 20 | 3.141623806667839 |
| 30 | 3.1416019000322795 |
| 40 | 3.141596556793832 |
| 50 | 3.1415946525910106 |
| 100 | 3.1415929035585526 |
We are as good as Aryabhata after 20 terms. In 40 terms, we are ahead of him by one digit. By term 100, we have one more digit of accuracy. Good, but we need to do better.
We note that the स्थौल्य is now
\( \begin{align*}
E_p &= \frac{1}{a_{p-1}} + \frac{1}{a_{p}} – \frac{1}{p} \\
&= \frac{1}{2(p-1)} + \frac{1}{2(p+1)} – \frac{1}{p} \\
&= \frac{1}{p^3-p}
\end{align*}\)
We note that a positive स्थौल्य indicates an over-correction, which we can verify from the table. It goes down as the third power of p, indicating that we now converge much faster with this correction term, as we can see from the table.
Second order correction
To improve the correction term, we therefore have to increase its denominator. Can we increase it by a constant 1? That will give us
\( \begin{align*}E_p &= \frac{1}{a_{p-1}} + \frac{1}{a_{p}} – \frac{1}{p} \
&= \frac{1}{2p-1} + \frac{1}{2p+3} – \frac{1}{p} \
&= \frac{-2p+3}{4p^3+4p^2-3p}
\end{align*}\)
We have now ended up with a term proportional to p in the error numerator. Therefore, Yuktibhāṣā reasons that the denominator of the end-correction is better increased by a factor of 1/p rather than 1, and arrives at the optimal value as \( a_p = 2p+2 + \frac{4}{2p+2}\) That gives us a स्थौल्य of
\( \begin{align*}E_p &= \frac{1}{a_{p-1}} + \frac{1}{a_{p}} – \frac{1}{p} \
&= \frac{1}{2p-2+\frac{4}{2p-2}} + \frac{1}{2p+2+\frac{4}{2p+2}} – \frac{1}{p} \
&= \frac{-4}{p^5+4p}
\end{align*}\)
This is a good स्थौल्य – negative, implying a slight under-correction, and reducing as the fifth power of p, which means we converge fast with this correction term. Let's see if that really happens
| Terms | Result with second order correction | Accuracy (decimal places) |
| 10 | 3.141590242370799 | 5 |
| 20 | 3.141592576186889 | 7 |
| 30 | 3.1415926433443224 | 7 |
| 40 | 3.1415926511540886 | 8 |
| 50 | 3.1415926527909903 | 9 |
| 100 | 3.1415926535648024 | 11 |
We have nine decimal places of accuracy with 50 terms, and 11 with 100, which is quite impressive. So, our correction term corresponds to
\( a_p = 2p+2 + \frac{4}{2p+2}\)
The term itself will be proportional to
\( \frac{1}{a_p} = \frac{1}{2p+2 + \frac{4}{2p+2}} = \frac{\frac{p+1}{2}}{(p+1)^2+1}\)
(ie: तस्या ऊर्ध्वगता या समसङ्ख्या तद्दलं गुणोऽन्ते स्यात् || तद्वर्गो रूपयुतो हारो …). Adding this correction term, our series now becomes:
\( C \approx 4.d – \frac{4.d}{3} + \frac{4.d}{5} – \frac{4.d}{7} + … (-1)^{\frac{p-1}{2}}.\frac{4.d}{p} + (-1)^{\frac{p+1}{2}}.\frac{4.d.\frac{p+1}{2}}{(p+1)^2+1} \)
The correction term has the opposite sign to the last term in the series – we add if we have subtracted, and vice versa (स्वमृणे कृते धने क्षेप एव करणीयः).
This form, as the table above attests, is good enough for calculation (लभ्धः परिधिः सूक्ष्मो), already to a degree of accuracy that was not known anywhere in Madhava's time. Comparing this table with the first, we see that it is a lot more useful to compute this correction term than to compute many many more terms of the series (बहुकृत्वो हरणतोऽतिसूक्ष्मश्च).
By improving the end correction term further, we can do even better than this, but we will leave that for a later article.
