BigInteger Square Root in F#

• Download SqrRootFSharp - 43.77 KB
• Download SquareRoot - 20.01 KB

Introduction

I have been dabbling around with the BigInteger type in F# and noticed that .NET does not provide an integer square root function for this type. So I have built a square root function that takes a BigInteger and returns a BigInteger. In F# parlance, it’s BigInteger -> BigInteger. One of my goals in this code is to only use the BigInteger type. The reason for this is that if you use some other type, say float, then you will restrict the size of number for which the function will work.

I also had my eye on performance since BigInteger math is notoriously slow. I use the Newton-Raphson method to converge to the result (no surprise there). The real performance gain is realized in how I calculate the first guess to seed the Newton-Raphson method.

Like many of you, I am new to F#; so any constructive suggestion to improve the code or to write in a more functional manner is welcome.

Background

The square root function is a frequently used function in all kinds of mathematical calculations, from finding the hypotenuse in the Pythagorean theorem to finding the distance between two points on a grid. It has many uses. Why it was left out of the library is anyone’s guess. Perhaps it will eventually be added to the .NET library.

Here is the definition of an integer square root: Given a positive integer S, return a positive integer R that is the greatest integer less than or equal to the square root of S.

This essentially says that if S is not a perfect square, then the answer is truncated down to the closest integer. If S is a perfect square, then you get the exact answer. For example, the integer square root of 25 is 5 and the integer square root of 26 is also 5.

Using the Code

We will start by building the Newton-Raphson method in F#. This is the standard function that is widely used by many calculators. I will not be explaining the Newton-Raphson method for finding the root; there are plenty of other sites where it is explained well. This method has a complexity on the order of O(Log n), that is log base 10; whereas the bisection method of finding the root has a complexity on the order of log base 2. Here is the main part of the code:

let bigintSqrt(bigNum : BigInteger) =
    let rec converge prevGuess =
        let nextGuess = (bigNum / prevGuess + prevGuess) >>> 1
        match BigInteger.Abs (prevGuess - nextGuess) with
        | x when x < 2I -> nextGuess
        | _ -> converge nextGuess
    if bigNum.Sign = -1 then
        failwith "Square root of a negative number is not defined."
    else
        let root = converge (sqrtRoughGuess bigNum)
        if root * root > bigNum then
            root - 1I
        else
            root		

I use the right shift operator (>>>) instead of dividing by 2 for a little better execution performance. The reason I added the last if statement (in the last four lines) is because it is doing an integer square root. So the answer must be less than or equal to the decimal square root. Since this loop will converge if it is within one, I may only need to change it by one. Also note that the internal function named converge is tail recursive so stack overflow is not a concern. (No other processing occurs after the recursive call.)

What’s left to do now is create the sqrtRoughGuess function used to seed the square root function. Since we are dealing with the BigInteger type, the input number could be hundreds or thousands of digits long, or more. Having a first guess of 1 would be a terrible guess in this case. Or having a first guess of our original number S or S/2 would be too high. Although any number greater than zero will work, the closer the first guess is to the eventual answer, then the function will do fewer iterations in the converge internal function.

In order to get a good first guess at the square root, I used a simple logarithmic identity in the following formula.

Sqrt(x) = b^{(1/2 * log(x))} , and b is the base of the logarithm in this formula. I’ll use b = 10 for the base in this formula to make things simple.

Here’s what I came up to get a close first guess.

let sqrtRoughGuess (num : BigInteger) =
    let log10x = log10 (num + 1I)
    let halfLog = (log10x + 1I) >>> 1  
    (power10 halfLog) 

I'm not showing the code in this article for the power10 function. But you can download the code from the link above to examine it. The function take 10 and raises it to the power of the given parameter. The power10 function also has a signature of BigInteger -> BigInteger.

Now, in order meet the criteria to use only the BigInteger type, I will create a log base 10 function. That is because the static Log10 method on the BigInteger type returns a float and to meet the goal to only use the BigInteger type.

let ten = 10I  // 10^1 
let tenBillion = 10000000000I // 10^10
let googol = 100000000000000000000000000000000000000000000000000000000000
00000000000000000000000000000000000000000I // 10^100
let googolToTenth = googol * googol * googol * googol * googol * 
googol * googol * googol * googol * googol // 10^10000

//*************************************************************************************
let private log10 theNumber =
    if theNumber <= 0I then
        failwith "Logarithm of a non-positive number is not defined."
    // this inner functions counts the number of times divisor will divide num
    let rec divisions count (num : BigInteger) divisor =
        match num with
        | x when x < divisor -> count, num
        | _ -> divisions (count + 1I) (num / divisor) divisor
    let thousandsCount, rem1 = divisions 0I theNumber googolToTenth
    let hundredsCount, rem2 = divisions 0I rem1 googol
    let tensCount, rem3 = divisions 0I rem2 tenBillion
    1000I * thousandsCount + 100I * hundredsCount + ten * tensCount + 
		fst(divisions 0I rem3 10I)

The log10 function, since it is an integer logarithm function, reverts to counting the number of times the integer parameter is divisible by 10. With this function complete, we have our square root function that only uses functions that take the BigInteger type. So the only bounds to how large the input integer can be is the capability of your computer.

Points of Interest

F# is an interesting language for someone who is an object oriented programmer and not used to the functional style. Recursion is the preferred form of looping rather than for / while loops. And tail recursion is necessary to prevent stack overflow errors. In the recursive functions in this article, I used accumulators in the list of parameters to avert other processing after the recursive call. You could also use continuations to ensure tail recursion, but I will leave that to the reader to research.

I also like the fact that I can mix the .NET languages in one solution. The GUI for this WinForms application is written in C#.

History

• 10 October, 2011 - First release
• 19 October, 2011 - Included executable download