Forums: » Register « |  Free Tools |  User CP |  Games |  Calendar |  Members |  FAQs |  Sitemap |  Support |
 User Name: Password: Remember me

New Free Tools on Dev Shed!
We're Excited to announce that Dev Shed now has 70 free tools on the site. To learn more, click here!

 Dev Shed Forums Sponsor:
#1
August 8th, 2011, 04:07 AM
 leszek31417
Contributing User

Join Date: Jul 2011
Posts: 313
Time spent in forums: 1 Week 2 Days 18 h 56 m
Reputation Power: 0
Integer square root algorithm

Hi,

i'm searching for efficient algorithm which computes integer square root of 64 bit (say, signed) integer number (long long).

That is: ISQ(n) = max k (integer) such that k*k <= n.

My approach is stupid.
I first compute k as integer part of sqrt(n) function (from <math.h>) and then try to increase or decrease k in a loop.

Thanks in advance for any suggestions.

Last edited by leszek31417 : August 8th, 2011 at 04:15 AM.

#2
August 8th, 2011, 05:46 AM
 salem
Contributed User

Join Date: Jun 2005
Posts: 4,264
Time spent in forums: 2 Months 4 Weeks 1 Day 17 h 18 m 9 sec
Reputation Power: 1827
__________________
If you dance barefoot on the broken glass of undefined behaviour, you've got to expect the occasional cut.
If at first you don't succeed, try writing your phone number on the exam paper

#3
August 8th, 2011, 08:19 AM
 leszek31417
Contributing User

Join Date: Jul 2011
Posts: 313
Time spent in forums: 1 Week 2 Days 18 h 56 m
Reputation Power: 0
Thank you.
Especially the second citation is something that i should read.

#4
August 9th, 2011, 12:10 PM
 mah\$us
Contributing User

Join Date: Feb 2009
Posts: 188
Time spent in forums: 3 Days 5 h 46 m 15 sec
Reputation Power: 49
This little Python function is not optimal, but I think is nonetheless pretty fast:

Code:
```def isqrt(n):
# construct an initial guess value
d = len(str(n))	# count of decimal digits
x = '3'
d = d / 2
while d > 0:
d = d - 1
x = x + '1'
x = int(x)

# perform Newton-Raphson iteration
prev = 0
while x != prev:
prev = x
q = n / x
x = (x + q) / 2

return x```
Python doesn't have an integer square root, or at least didn't back when I wrote this function.

The first block makes a crude approximation to the square root; the second repeats Newton-Raphson until it converges. For extra credit, prove that the while loop will always terminate (or, prove that it can fail to terminate)

#5
August 10th, 2011, 04:57 AM
 leszek31417
Contributing User

Join Date: Jul 2011
Posts: 313
Time spent in forums: 1 Week 2 Days 18 h 56 m
Reputation Power: 0
Thank you again !

I'll play with this algorithm.

How many steps the second loop (Newton-Raphson iteration) needs?

#6
August 11th, 2011, 03:22 AM
 mah\$us
Contributing User

Join Date: Feb 2009
Posts: 188
Time spent in forums: 3 Days 5 h 46 m 15 sec
Reputation Power: 49
Quote:
 How many steps the second loop (Newton-Raphson iteration) needs?
Newton-Raphson for square root is a wonderful algorithm, ranking in my mind with Euclid's. It is both very simple, and extremely efficient. When my crude technique is used for the initial "guess" value, the maximum number of iterations required to reach the nearest integer square root of n is well approximated by log2(log2(n)). [log2 means "log base 2", or roughly the number of bits.]

Because the program I posted uses one more iteration to detect convergence, a good estimate is 1 + log2(log2(n)). Note well that in some cases, the number of iterations will be much lower than this worst-case figure.

Numbers up to 10,000 bits require no more than about 13 iterations; up to 100,000 bits, no more than about 16 iterations.

As I hinted in the previous post, the Python script I gave doesn't always terminate: in integer math, the iteration can stably alternate between two values, in which case the test for same value as before will never evaluate true. Here's a better version:

Code:
```def isqrt(n):
# restrict to safe arguments
if n < 1: return 0

# construct an initial guess value
d = len(str(n))	# count of decimal digits
x = '1'
d = d / 2
while d > 0:
d = d - 1
x = x + '0'
x = int(x)

# perform Newton-Raphson iteration
prev = 0
prv2 = 0
while x != prev and x != prv2:
prv2 = prev
prev = x
q = n / x
x = (x + q) / 2

return x```
It can be improved by computing a more refined initial value for Newton-Raphson, which could reduce the worst case by probably one or two iterations.

#7
August 12th, 2011, 02:53 AM
 leszek31417
Contributing User

Join Date: Jul 2011
Posts: 313
Time spent in forums: 1 Week 2 Days 18 h 56 m
Reputation Power: 0
Thank you very much Mr. Mah\$us !

I played with this algorithm a little.

Quote:
 Originally Posted by mah\$us ... is well approximated by log2(log2(n)) ...

Where i can find the proofs of the above?

Here is my C version (without good initial guess).

PHP Code:
```  int isqrt ( int n ){ int p,r;        r = n; // initial guess    do{               p = r;       r = ( r + n/r )/2; // integer division!           }while( r < p );        return p; }  ```

#8
August 12th, 2011, 01:46 PM
 mah\$us
Contributing User

Join Date: Feb 2009
Posts: 188
Time spent in forums: 3 Days 5 h 46 m 15 sec
Reputation Power: 49
1. Though the original post specified 64 bits, and integers this size are easily handled in modern C implementations, computation related to cryptography sometimes requires taking square roots of much bigger numbers (more than 1000 bits).

Anyone who wants to write a square-root program in C for large integers can use a large-integer library (there are several available). Python natively handles large integers, and so is very convenient for this kind of computation.

2. leszek, the loop termination in the C program is valid, if and only if the successive approximations decrease monotonically. I suggest that you prove (or disprove) this monotonicity.

3. The proof for the number of iterations is not difficult, and rather fun, so I'll leave it as an exercise. Hint: if y is the true square root, define the most recent approximation as y (1 + epsilon), where epsilon = 2^-n. Speaking a little loosely, we can say that this approximation is "good to n bits" (that is, the first n bits are correct). Using the iteration formula to compute the next approximation -- that approximation is good to how many bits? Once you know this, it is not difficult to compute the number of iterations required.

#9
August 13th, 2011, 03:43 AM
 leszek31417
Contributing User

Join Date: Jul 2011
Posts: 313
Time spent in forums: 1 Week 2 Days 18 h 56 m
Reputation Power: 0
second version

PHP Code:
``` // second version int isqrt ( int n ){ int p,r;    p = 1; while ( p <= n/p ) p = 2*p;        r = p; // initial guess: isqrt(n) < r        do{          p = r;       r = ( r + n/r )/2;           }while( r<p );        return p; }  ```

It works very fast! Now i believe in log2(log2(n)) convergence!

I am now working on proof of correctness.

#10
August 16th, 2011, 03:31 AM
 leszek31417
Contributing User

Join Date: Jul 2011
Posts: 313
Time spent in forums: 1 Week 2 Days 18 h 56 m
Reputation Power: 0
Let k be the integer root of n (i.e., k is the number we are looking for).
Let r' = ( r + n/r )/2 (integer divisions and r>0).

Lemma. r'<r if and only if r>k.
Proof. Easy.

So the sequence of iterations is strictly decreasing, provided it stays ABOVE k.

This proves the termination of algorithm: if we start with initial value > k,
then there must be a step such that r' <= k < r .
The next iteration value (say r'') is greater than or equal to r' (from the "only if" part),
so the loop terminates.

Claim. r' = k.
Proof. A little bit more involved. As main tool use the well known fact
that geometric mean <= arithmetic mean.

Mr. Mah\$us, thank you very much for bringing to my attention this beautiful algorithm !

Last edited by leszek31417 : August 17th, 2011 at 02:32 AM. Reason: english language ...

 Viewing: Dev Shed Forums > System Administration > Security and Cryptography > Integer square root algorithm

## Developer Shed Advertisers and Affiliates

 Thread Tools Search this Thread Search this Thread: Advanced Search Display Modes Rate This Thread Linear Mode Rate This Thread: 5 : Excellent 4 : Good 3 : Average 2 : Bad 1 : Terrible

 Posting Rules You may not post new threads You may not post replies You may not post attachments You may not edit your posts vB code is On Smilies are On [IMG] code is On HTML code is Off
 View Your Warnings | New Posts | Latest News | Latest Threads | Shoutbox Forum Jump Please select one User Control Panel Private Messages Subscriptions Who's Online Search Forums Forums Home -------------------- Programming Languages    PHP Development        PHP FAQs and Stickies    Perl Programming        Perl FAQs and Stickies    C Programming        C Programming FAQs and Stickies    Java Help        Java FAQs    Python Programming        Python Programming FAQs    Ruby Programming        Ruby Programming FAQs    Game Development        Game Development FAQs Programming Languages - More    ASP Programming        ASP Programming FAQs    .Net Development        .Net Development FAQs    Visual Basic Programming        Visual Basic Programming FAQs    Software Design        Software Design FAQs    ColdFusion Development        ColdFusion Development FAQs    Delphi Programming        Delphi Programming FAQs    Regex Programming        Regex Programming FAQs    XML Programming        XML Programming FAQs    Other Programming Languages        Other Programming Languages FAQs Web Design    HTML Programming        HTML Programming FAQs    JavaScript Development        JavaScript Development FAQs    CSS Help        CSS Help FAQs    Flash Help        Flash Help FAQs    Photoshop Help        Photoshop Help FAQs    Web Design Help        Web Design Help FAQs    Website Critiques        Website Critiques FAQs    Search Engine Optimization        Search Engine Optimization FAQs Mobile Programming    Mobile Programming        Mobile Programming FAQs    iPhone SDK Development        iPhone SDK Development FAQs    Android Development        Android Development FAQs    BlackBerry Development        BlackBerry Development FAQs Web Site Management    Business Help        Business Help FAQs    Development Software        Development Software FAQs    Scripts        Scripts FAQs Databases    Database Management        Database Management FAQs    DB2 Development        DB2 Development FAQs    MySQL Help        MySQL Help FAQs    PostgreSQL Help        PostgreSQL Help FAQs    Firebird SQL Development        Firebird SQL Development FAQs    MS SQL Development        MS SQL Development FAQs    Oracle Development        Oracle Development FAQs    LDAP Programming        LDAP Programming FAQs System Administration    Mail Server Help        Mail Server Help FAQs    Apache Development        Apache Development FAQs    Security and Cryptography        Security and Cryptography FAQs    Antivirus Protection        Antivirus Protection FAQs    DNS        DNS FAQs    IIS        IIS FAQs    Networking Help        Networking Help FAQs    FTP Help        FTP Help FAQs Operating Systems    BSD Help        BSD Help FAQs    Linux Help        Linux Help FAQs    UNIX Help        UNIX Help FAQs    Windows Help        Windows Help FAQs    Mac Help        Mac Help FAQs Web Hosting    Web Hosting        Web Hosting FAQs    Free Web Hosting        Free Web Hosting FAQs    Web Hosting Requests        Web Hosting Requests FAQs    Web Hosting Offers        Web Hosting Offers FAQs Computer Hardware    Computer Hardware    CPUs        CPUs FAQs    Cooling        Cooling FAQs    Embedded Programming        Embedded Programming FAQs    Motherboards        Motherboards FAQs    Multimedia Hardware        Multimedia Hardware FAQs Other    Dev Shed Lounge        Dev Shed Lounge FAQs    Development Articles        Development Articles FAQs    Beginner Programming        Beginner Programming FAQs    Hire A Programmer        Hire A Programmer FAQs    Project Help Wanted        Project Help Wanted FAQs Latest News Updated Hourly    Technology News    Business News    Science News Forum Information    Forum Rules/Guidelines        Forum Rules/Guidelines FAQs    Forum Announcements        Forum Announcements FAQs    Dev Shed Gaming Center        Go to the Dev Shed Battle Arena        Go to the Dev Shed Arcade Games        Go to the Legend of the Green Dragon    Suggestions & Feedback        Suggestions & Feedback FAQs

 Forums: » Register « |  Free Tools |  User CP |  Games |  Calendar |  Members |  FAQs |  Sitemap |  Support |