Computing a random root of a number

To compute a random root of a number is no trivial task. In a computer this is normally done in an iterative function.

The equation


Roots

is switched to


Roots

and the root for this equation is searched. This root finding is done in the similar manner as is done in a Newton/Maehly algorithm (see Newton/Maehly )


We start at a particular point and run along the function till it crosses the cero line. In the computer program this is done by starting at point P1 and laying a tangent onto the graph. This tangent should cross the cero line somewhere at x2 and there where it crosses it is already a bit closer to the cero crossing of the graph than where we started. From x2 we go to P2 and do the same again: Lay another tangent onto the graph and go to its cero crossing x3 and do the same again till we are close enough at the cero crossing. Close enough means the difference between one such iteration and the last is small enough.


Roots

In mathematics that means with


Roots

we get


Roots

and so on… same as in the Newton/Maehly algorithm. Only the function is much simpler.

As


Roots

we can set


Roots

With this formulation we start with x = radiant and run down the function. For the square root of 9 that would look like



double radiant = 9;
double exp = 2;
double x_old = radiant;
double x_new = x_old;
bool run = true;
 
while (run)
{
     

x_new = x_old;

for(int i = 1; i < exp-1; i++)

x_new = x_new * x_old;

x = x_old + (radiant - x_new * x_old) / exp / x_new;

run = Math.Abs(x_old - x) > 1E-8;

x_old = x;

}



And that gives x = 3.00 at the end of the loop.


Now there are some special cases that need to be handled:

This loop only works if the exponent is bigger than 1. If the exponent is 1, the result of x1 will be x. If exp is negative, we have 1/ x-exp or finally if exp = 0, the result is just 1.

All these special cases implemented in a function Root(double radiant, double exp) would look like:



public static double Root(double radiant, int exp)
{
double x = 0;
double x_old = radiant;
double x_new = x_old;
bool run = true;
int count = 0;
if (exp > 1.0)
{
while (run && (count < 1000))
{
x_new = x_old;
for(int i = 1; i < exp-1; i++)
x_new = x_new * x_old;
x = x_old + (radiant - x_new * x_old) / exp / x_new;
run = Math.Abs(x_old - x) > 1E-8;
x_old = x;
count++;
}
if (count >= 1000)
x = 0;
}
else
{
if (exp == 1)
x = radiant;
else
{
if (exp < 0)
{
x = Root(radiant, -exp);
if (x != 0)
x = 1 / x;
}
else // exp is 0
x = 1;
}
}
return x;
}



It returns the result if it is valid and 0 if it is not valid.

Like this a computer calculates a random root of a random number.

I use this function to compute the square root of 3 with 580 decimal places in Computing Pi :-)