Chaos in the iterative Hindu square root method of the gaṇaka-rāja

For Hindus big numbers always mattered and our mathematics is quite reflection of this fascination. Since the earliest times, Hindus devised various methods to obtain square roots of numbers, especially approximations of irrational roots correct to multiple decimal places. The earliest of these methods involving a series of terms is seen encoded in the altars for the Soma rituals specified in the saṃhitā-s of the Yajurveda and explicitly spelled out in their the śulbasūtra-s. Indeed, we have evidence that development of these methods continued in the Yajurvaidika tradition as indicated by Rāma dīkṣita’s commentary on Kātyāyana where he provides a tradition regarding a further term to the approximation to get \sqrt{2} correct to 7 decimal places. A similar improvement was likely used in the procedure preserved by Sundararāja dīkṣita in the Āpastamba tradition for an approximate squaring of the circle based on \sqrt{2}.

By the last few centuries before the common era the Hindus had already discovered a method similar to what is today known in the west as the first term Newton-Raphson approximation. We also see the exact algorithm for both square roots and cube roots of ācārya Āryabhaṭa further explained for the lay by Bhāskara-I. But the high point of the Hindu tradition of iterative methods is seen in the text of the brāḥmaṇa Chajjaka-putra gaṇaka-rāja probably from Mārtikāvati (unfortunately named Bakshali manuscript: BM), which gives a glimpse of just what Hindu knowledge has been lost over the ages. While this method was misunderstood by the earlier white indological translator of the BM, the sophistication of the gaṇaka-rāja’s method has only more recently become clear. This has been explained and commented upon in detail by the computer scientists Bailey and Borwein in their excellent work on the same. We shall here comment upon an interesting aspect we discovered of the functions involved in the method .

While the method has already been discussed in detail by Bailey and Borwein, we shall go over it here for introducing the system. In order the find the square root of a number q the BM suggests the following procedure:
Take some starting number: x_n = x_0

Then y_n \approx \sqrt{q}. Now if we take x_n=y_n and iterate the above procedure we get increasingly accurate approximations of \sqrt{q}.

As a example let us take q = 5 and x_0=0.1. Then we have the following:
1)\; 12.6248003992015985\\ 2) \; 3.6392111847990769\\ 3) \; 2.2506636482615887\\ 4) \; 2.2360679780006203\\ 5) \; 2.2360679774997898
Thus, in iteration 3 the value of \sqrt{5} correct to 1 decimal place, in iteration 4 it is correct to 8 decimal places and in iteration 5 it is correct to at least 16 decimal places, in line with the Hindu love for big numbers.

Now if we instead take x_0=2 because we know that \sqrt{5} should lie somewhere in the vicinity of 2 then we get:
1) \; 2.2361111111111112\\ 2) \; 2.2360679774997898
Thus, with this close value right in the first iteration we get it correct to 3 decimal places and in the second to at least 16 decimal places! As Bailey and Borwein had shown it quartically converges on the square root. Now if we take a negative number for x_0 it then converges similarly to -\sqrt{q}.

hindu_figure1Figure 1

Now consider the following alternative procedure where instead of plugging x_n=y_n we plug x_n= x_{n+1} and thus generate for each iteration (x_{n+1},y_n). On plotting the map of (x_{n+1},y_n) we see the points fall on an interesting curve (Figure 1). This curve has two boat-shaped branches which are respectively tangential to the lines y= \pm \sqrt{q}. The region of tangency is peculiar in that the curve lingers in the proximity of y= \pm \sqrt{q} over a wide x-interval.


Figure 2

The actual map of the points obtained by the above procedure displays an interesting feature: they are spread all over the two branches of the curve above but fall most frequently in the vicinity of the two root lines. They notably decrease in frequency as one moves away from those lines but we do get to see extreme points far away from the two root lines. Thus, \pm \sqrt{q} serve as the peaks (Figure 2) for the distribution of y_n with a clear decline for greater and lesser allowed values respectively. However, the tails of their distribution are prominent enough that we seen multiple extreme values. The median value for the negative side of the distribution is \approx -2.3 and for the positive side is \approx 2.3, illustrating the dominance of the values close to \pm \sqrt{q}. In line with this, for a large enough number of iterations with a given x_0 the overall median value of y_n comes out as \pm \sqrt{q} . However, below is an examination of the extreme values reached by y_n for a run initiated with q=5, x_0=0.1 for 2000 iterations:
Minimum: -438.98149
Maximum: 133.19996
This shows that y_n explores values over 50-100 times the median values in course of the iterations.

To understand this map better let us look at it geometrically (Figure 3). The two expressions that are deployed successively by Chajjaka-putra to get the square root represent the below functions:
f(t)= \dfrac{q-t^2}{2t}\\ g\left(t\right)=t+f\left(t\right)-\dfrac{f\left(t\right)^2}{2\left(t+f\left(t\right)\right)}\\[10pt] \therefore g\left(t\right)=\dfrac{q^2+6qt^2+t^4}{4qt+4t^3}


Figure 3

We see that f(t) is a hyperbola with the y-axis as one of its asymptotes. g(t) is a quartic curve, which has y= \pm \sqrt{q} as the as its tangents with the points of tangency being (\sqrt{q},\sqrt{q}) and (-\sqrt{q},-\sqrt{q}). This curve has a very “flat” type of tangency, i.e. it lingers in the proximity of y= \pm \sqrt{q} over an extended x-range. This is the secret of the gaṇaka-rāja’s method firmly “pulling” things to the vicinity of required square root. Thus, the parametric curve (f(t),g(t)) is the one on which the points of the above-described map based on gaṇaka-rāja’s two expressions lie (Figure 1, 3). This curve as noted from the above map has y= \pm \sqrt{q} as its tangents with the point of tangency at (0,0), but like g(t) it lingers over a wide x-range close to the point of tangency. This explains why the map tends to concentrate the points in the vicinity of y=\pm \sqrt{q}.

Now, if we look at the positions of actual points of the map on (f(t),g(t)) an interesting observation becomes apparent: while tending to cluster in the vicinity of y= \pm \sqrt{q}, successive points are not necessarily proximal to each other on the curve. Rather they jump about the curve in a chaotic fashion (Figure 3) either on the same branch or between the two branches. This becomes even rather apparent if we plot a run of y_n against iteration number n=1..2000 (Figure 4; The gaps in the plot are where the y_n jumped outside the range of \pm 50 which we set for good visualization). What is notable is that the rarer values in addition to appearing chaotically are also rather extreme: the sum of positive y_n for 2000 iterations with x_0=0.1 is 4578.39; of this just 12 values add up to 1005.48, which is \approx 22 percent of the total sum. Each of these 12 extreme values is over 20 times the median value of positive y_n. The picture is roughly symmetric for the negative values.


Figure 4

Importantly, this chaotic behavior of y_n is very sensitive to the initial values of x_n with which we start the map. This is dramatically illustrated by two points close to root of q=5, 2.2 and 2.21 (Figure 5): while for the first 9 iterations the evolution of the two initial values is the same in direction though not magnitude, iteration 10 and beyond they go completely out of synchrony in both direction and magnitude.


Figure 5

In conclusion this map provides an analogy to think about certain processes in nature and historical events. First, it provides a potential model for foraging behavior of variety of organisms. In this model the clustering around the root values represents what might be called the base-line or ordinary foraging and the extreme jumps represent the drastic forays away from their local patch to distant locales. Such behavior may be seen in animals among herbivores moving to new feeding grounds far from their usual feeding areas or certain carnivores like sharks seeking new hunting waters far from their their current zone. This kind of behavior is also seen in certain ciliates like Halteria in the microscopic realm. In the world of protein sequences we see a similar tendency to keep to a tightly constrained space of diversity under purifying selection within which there is a low-radius exploration under neutral drift. This is punctuated by huge saltations that result from strong positive selection for new functional niches. This might happen within a family of proteins or in the proteome of an organism with some proteins showing big saltations in sequence space.

The great mathematician Benoit Mandelbrot, the pioneer in the study of chaos, has brought home the importance of distributions with rare events with extreme values. The role of such events in systems like financial markets has been recently explained at length by Nassim Taleb. In this context, the above-described map provides an analogy for one type of historical evolution of systems. Even with same basic parameter ( \sqrt{q}), clearly predictable bulk statistics (e.g. median and range of most frequent values) and similar starting values we see: 1) clear differences in long term evolution with non-overlapping chaotic extreme events different in both magnitude and direction. 2) Extreme events that can disproportionately contribute to the total numerical measure of the events in the series. A historical system evolving under such a model shows us how with very similar starting material and bulk behavior we can have a great difference in actual events and outcomes. This might be similar to actually observed phenomena like the fall empires or the sudden extinction of long-lasting lineages. This is a theme which we might explore further with examples of some other such systems which have been studied for their chaos.

~ by mAnasa-taraMgiNI on October 21, 2016.

%d bloggers like this: