# What is a fast algorithm for finding the integer square root?

What is a fast algorithm for computing integer square roots on machines that doesn't support floating-point arithmetic?

I'm looking for a fast algorithm for computing the integer square root of an integer $$N \ge 0$$, i.e., the smallest integer $$m \ge 0$$ such that $$m^2 \le N < (m+1)^2$$.

The reason is that I'm operating on a virtual machine that doesn't support floating-point arithmetic (real-numbers), so algorithms like Newton's cannot be implemented, right?

The naive solution for finding the square root of $$N$$ is to check for $$m=0,1,\ldots$$ whether $$m^2 \le N < (m+1)^2$$. Is this the best algorithm?

• Apart from the integer version of Newton's algorithm, a binary search is also faster than the naive method (except for very small $N$, but then all methods that aren't insane are fast). Oct 12, 2017 at 18:10
• The optimal solution highly depends on the amount of RAM and machine commands available, and whether are you going to perform arbitrary precision operations or just use a fixed-size (which one) integer. Jun 24, 2020 at 10:19

I think the easiest method will Babylonians sqrt method and it works well with machine supporting only integer division.

def babylonian(n):
x = n
y = 1
while(x > y):
x = (x+y)//2
y = n//x
return x


This algorithm runs in approximately $$O(\log n)$$ time. More details can be found on StackOverflow.

A fast ($$O(\log n)$$) way to calculate the integer square root is to use a digit-by-digit algorithm in base2:

$$\text{isqrt(n)} = \begin{cases} n & \text{if n < 2} \\ 2\cdot\text{isqrt(n/4)} & \text{if (2\cdot\text{isqrt(n/4)} + 1)^2 > n} \\ 2\cdot\text{isqrt(n/4)+1} & \text{otherwise} \end{cases}$$

Making sure to calculate $$n/4$$ using bitshifts, and doing the above iteratively. An example for 16-bit integers can be found on Wikipedia.

Assuming integer log2 is available, you can get a pretty good starting estimate by rewriting $$\sqrt{x}$$ like this: $$\sqrt{x}=x^{\frac{1}{2}}=2^{\frac{\log_{2}{x}}{2}}$$ Integer $$\log_{2}$$ is often one or two cheap instructions. You just have to refine the result a couple times since everything is floored. This O(1) implementation in C++ finds the square root of every integer $$[0, 2^{32}-1]$$ in 15 seconds on my PC - about 3.5 ns per iteration.

#include <bit>
constexpr uint16_t Sqrt(const uint32_t x) {
// Avoid divide by zero
if (x < 2) {
return static_cast<uint16_t>(x);
}
// This code is based on the fact that
// sqrt(x) == x^1/2 == 2^(log2(x)/2)
// Unfortunately it's a little more tricky
// when fast log2 is floored.
const uint32_t log2x     = std::bit_width(x) - 1;
const uint32_t log2y     = log2x / 2u;
uint32_t       y         = 1 << log2y;
uint32_t       y_squared = 1 << (2u * log2y);
int32_t        sqr_diff  = x - y_squared;
// Perform lerp between powers of four
y += (sqr_diff / 3u) >> log2y;
// The estimate is probably too low, refine it upward
y_squared = y * y;
sqr_diff  = x - y_squared;
y += sqr_diff / (2 * y);
// The estimate may be too high. If so, refine it downward
y_squared = y * y;
sqr_diff  = x - y_squared;
if (sqr_diff >= 0) {
return static_cast<uint16_t>(y);
}
y -= (-sqr_diff / (2 * y)) + 1;
// The estimate may still be 1 too high
y_squared = y * y;
sqr_diff  = x - y_squared;
if (sqr_diff < 0) {
--y;
}
return static_cast<uint16_t>(y);
}


The two variable integer divisions in the refinement step aren't cheap, but this was still faster than any of the looping solutions I came up with.

For $1 \leq N \leq 15,$ suggest table lookup. $1 \leq N \leq 3,$ output $1.$ For $4 \leq N \leq 8,$ output $2.$ For $9 \leq N \leq 15,$ output $3.$

Otherwise, as far as finding a good starting point, I suppose you can begin with $1,4,16,..., 4^k$ until $4^k \geq N > 4^{k-1}.$ Then let the starting point be $$x_0 = 2^k$$

Newton's method with integers and the floor function, but at the very end prevent loops: $$x_{j+1} = \left\lfloor \frac{x_j^2 + N}{2 x_j} \right\rfloor$$ As soon as $$|x_{j+1} - x_j| < 5,$$ switch to your one-at-a-time idea between those endpoints. Otherwise, an infinite loop is possible.

When I programmed this, I let Newton go until consecutive terms were very close together ($<5$). Then I set $m$ to the smaller one, and increased $m$ by one until the square exceeded $N.$ Then I decreased $m$ by one until the square is no larger than $N.$

• You have found the integer square root - namely $x_j$ - when $x_{j+1} \geqslant x_j$ with $j \geqslant 1$ (that is in case the initial guess $x_0$ was too small, then $x_1 > x_0$ doesn't imply that $x_0$ is the integer square root). Oct 12, 2017 at 18:06
• This is slow because it uses integer division.
– orlp
Oct 12, 2017 at 18:09
• @DanielFischer I did not put in everything. When I programmed this, I let Newton go until consecutive terms were very close together. Then I set $m$ to the smaller one, and increased $m$ by one until the square exceeded $N.$ Then I decreased $m$ by one until the square is no larger than $N.$ Oct 12, 2017 at 18:17
• I'm just saying that there is a quite simple stopping rule. Iterate until $x_{j+1} \geqslant x_j$, then $x_j$ is it (if $j > 0$). Oct 12, 2017 at 18:20
• @DanielFischer I see what you mean now; that seems better than what I did. At the time, I was dismayed to find an infinite loop I did not understand, until I started putting in cerr commands and so on. The thing I described was my cure, first thing to come to mind. Oct 12, 2017 at 18:23