Integer Square Roots

Well, I don't know the answer, but I might be able to help you with giving you a subquestion. Square root is nothing else then: X^(1/2). How many stepps (cycles) do you need to take to perform this operand, with the way a x86 processor do there math operations.

I know I'm very vague, but I really don't know the aswer. I hope this helped you a bit anyway.

Quote:

Original post by Xeile Well, I don't know the answer, but I might be able to help you with giving you a subquestion. Square root is nothing else then: X^(1/2). How many stepps (cycles) do you need to take to perform this operand, with the way a x86 processor do there math operations. I know I'm very vague, but I really don't know the aswer. I hope this helped you a bit anyway.

well, the answer to your question is very dependant on the algoritm you pick, and i dont know which to pick.

what i currently do is this:

traverse the bytes of the outcome in significance decending order. if setting this bit wont raise the square of the outcome above the input, do so. that, unrecognizably simplified and completely unrolled for all bits.

per bit its:

1 compare 1 shift 1 add1 or1 sub

so i suppose 16bits * 5 (max) = 65 one cycle operations isnt bad: atleast the compiler settings are not the problem. however, its the algorithm itself im concerned with.

WHOA.

who was the guy telling me 64 bit operations arent slower than 32 bit ones on 32 bit platforms?


this runs in the aforementioned 65 cycles.

only making these changes:

ADDS 100 CYCLES!

*me runs to the store to buy a 64 bit processor*

Quote:

Original post by Anonymous Poster wtf is wrong with int x = (int)sqrt(y); ?

nothing. its more of a personal venture into the world of fixed point math than anything 'usefull'.

Quote:

I doubt that sqrt() is going to be your bottleneck.

i doubt that im interested in the standard stfunoob response. this is going to be called millions of times a second. although speed isnt my main goal, within the constraint of fixed point math is still consider it important.

Well, you don't want to perform square roots millions of times per second. Integer or not. ;)

Square roots *are* slow. Not much can be done about it, although in many situations, you can avoid performing them in the first place, through a bit of trickery.

On a side note, there's generally not much saved by using fixed point rather than float on modern processors.

Hi,

if you haven't found this site already, it's a nice collection of squareroot-algos:

http://www.azillionmonkeys.com/qed/sqroot.html

Quote:

Original post by Spoonbender Well, you don't want to perform square roots millions of times per second. Integer or not. ;)

true, but i dont have a choice.

Quote:

Square roots *are* slow. Not much can be done about it, although in many situations, you can avoid performing them in the first place, through a bit of trickery.

not here.

Quote:

On a side note, there's generally not much saved by using fixed point rather than float on modern processors.

i dont expect any savings, except for a bit of memory.

Quote:

Original post by h3idi Hi, if you haven't found this site already, it's a nice collection of squareroot-algos: http://www.azillionmonkeys.com/qed/sqroot.html

thanks, ill check it out.

Quote:

Original post by Eelco Quote: Original post by h3idi Hi, if you haven't found this site already, it's a nice collection of squareroot-algos: http://www.azillionmonkeys.com/qed/sqroot.html thanks, ill check it out.

ah yes ive already seen that. usefull page, too bad im not much ofan ASM guru, nor do they give a clear view on what is the fastest solution. ah well, i guess i can live with 200 cycles for a 32bit->32bit sqrt until i find something better.

If this will be runned millions of times per second and if the values of your integers are not gigantic, you can just make a table of the square roots, and look them up in constant time when you need them. So sqroot[64]=8 ect.

You could try a lut with the 16 bit square numbers (theres only 255 of them, so nothing to worry about). Or you could go and make a 64K lut (which would be very small, so no need to worry much about memory), which contains the square root of each 16 bit number.

You go and you can do a neat pointer trick to extract both parts of the number.

You go and you set up a pointer to an int16, and you set it to point to the number.

You then acess it as pointer[0] = The whole part pointer[1] is the decimal part.

Neat, huh?

Theres another way to go and find you the square root of your number.

You go and try this:

sqrt = ((n + 1) + ((n + n))) / (3 + n)

Taken from mathworld Here This should prove adiquite, seeing as your only looking for an integer square root.

Its only 4 additions (one paired so they can be done simultaniously), and a division. There might be a better algorithm somewhere.... You could try this when memory ops are expensive (like in desktop CPU's where a memory load takes more time then a few hundred clock cycles).

You can then do a liniar aproximation of the square root, of the integer part.

Thats basically
squ = (sqrt * fn) + ((Sqrt + 1) * (1 - fn))

Assuming fn is in the range from 0 to 1. Squ is then your answer.

Any more questions's?

I'll try and see if i can find a better answer somewhere...

ok, most of that was for finding the root of a 16.16 number...

One that you could use, would probably be a two stage process.
First one, using somehting like:

register int tmp;register int tmmp;register int tmmmp;tmp = 10*n; //Or you could use the bitshifts if you wanttmmp = n * n;tmmmp = tmmp + tmmp;Fstagesqrt = (1 + tmp + tmmmp + tmmmp + tmmmp) / (5 + tmp + tmmp);

Notice how that it uses mostly additions (which are verry easily parellised, and are very fast), and only two multiplications? (one if you use the two (parellised) bitshifts to make tmp).

Now this should be good enough for most people. (It should have an amazingly low error rate, because its the equivilent of three iterations of the algorithm.)

If that isn't good enough then you'll need to call in the big guns.

(1 + 28 * n + 70 * (n * n) + 28 * (n * n * n) + (n * n * n * n)) / ((8 + 8*n) * (1 + 6*n + (n * n))

Thats a bit harder to optimise, but it is very good. you could also and use a simpler algorithm which uses a guess (which could be the result of the other two square root algorithms) to converge onto an answer.

The algorithms simple:

Its simply ((G * G) + N) / (G + G)
That will work two...

Use whatever you need.

G's your guess, N's your number. And the result is an aproximate square root for the number.

From,
Nice coder

[Edited by - Nice Coder on January 28, 2005 2:08:01 AM]

thank you nice coder, for taking the time to write all that.

i am however only interested in the exact outcome, and a lookup table doesnt sound too appealing i must say. i dont want to trash my cache.

i tried extending my 32->16 bit squareroot algo to 32->16.16

however, it doesnt work and i dont get for the life of me why. what i simply do is shift the 32 bit input into the high bits of a 64bit. then apply the exact same but twice as long algo on this number, and interpret the 32 bit outcome as 16.16. however, the results are complete gibberish.

Eelco- Your wealcome, i could probably do a bit more, if you need to. (i'm not very good at fixed point, sorry...)

Yeha, fixed point square roots are very..... difficult.

Check the binary representation of the number, then of the expected number.
You could be suffering from a bp misplacement.

LUT's are good if your running on something that has slow processing vs. Speed. With modern cpu's its actually much faster (by a lot) just to compute the square root the normal way.

I would say, go and do your normal binary square root funtion, but go and rightshift the output number (because shifting the input number isn't too good for this), to go and get the 16.16 notation.

From,
Nice coder

Quote:

Original post by Nice Coder Check the binary representation of the number, then of the expected number. You could be suffering from a bp misplacement.

thats what i did, but sofar it hasnt given me any clues. its as if longs arent behaving like they should, but thats probably not so.

Quote:

I would say, go and do your normal binary square root funtion, but go and rightshift the output number (because shifting the input number isn't too good for this), to go and get the 16.16 notation. From, Nice coder

i have to shift the input number.

if you look at the algorithm i posted above, youll see why. it takes an n-bit input, and yields and n/2-bits output, which makes sense. so in order to end up with 16.16/32 (the same thing, its a matter of interpretation), ill need to input a 64 bit number.

but interpretations aside, it just seems like i cant scale my 32->16 algo to 64->32, even though it should be straightforward. well im going radical: just printing every intermediate values bitpattern to the screen, that should show me the problem :).

Ahem, are you sure, that for the 64 bit number your feeding it (not the 32 bit one with the bp displacement), that your not getting the correct results?

DENC

Quote:

Original post by Eelco +-65 cycles to take the 16 bit root of a 32 bit number. (set the highest possible value btw. for low values its something like 50).

Sure it's quite a pity.

Quote:

if i want the root in 16.16 fixed format, which i think i will, thats ofcource going to cost me twice as much time.

- Unless you adapt the algo. But in theory yes, since it makes 32 bits. And usually the order of the series progress linearilly with the number of bits you need.

- Also are you sure you really need a full 32 bits f precision ? For rendering purposes I doubt, since the output is usually 8 bits/color component.

- Using SIMD technologies a sqrt would cost you a few portions of cycles. Not 200 !

=> Thus, I'd not pay more than 1 cycle if it was for typical rendering.

Quote:

Is this bad? I have the feeling it can be better. i dont want to resort to any float ops though. its quite time critical, because it will be used for normalization in a raytracer on a per-ray basis.

".. dont want to resort to any float ops though => its quite time critical ..."

I hope that you did not mean an implication here. Because FPU is damn much speedier than integer arithmetics unless you consider SIMD instructions.

There is far enough room for improvement if it's really "time critical". So be certain that your prerequisites are really valid and motivated.

Else you can also boost it using LUTs then NR. In fact it's how fast 3DNow or SSE rsqrt (thus sqrt=x*rsqrt) work internally. But the nature of floating point encoding (exponent and mantissa) makes the LUT approach very astute and efficient.

So once again, are you certain that you don't want any FPU or SIMD in your routine ? I don't mean necessarilly the sqrt of "math.h", the "noob" implementation.

http://astronomy.swin.edu.au/~pbourke/analysis/sqrt/index.html

i was hoping you would reply charles :)

Quote:

Original post by Charles B Quote: Original post by Eelco +-65 cycles to take the 16 bit root of a 32 bit number. (set the highest possible value btw. for low values its something like 50). Sure it's quite a pity.

what would you consider good performance for an exact integer root?

Quote:

Quote: if i want the root in 16.16 fixed format, which i think i will, thats ofcource going to cost me twice as much time. - Unless you adapt the algo. But in theory yes, since it makes 32 bits. And usually the order of the series progress linearilly with the number of bits you need. - Also are you sure you really need a full 32 bits f precision ? For rendering purposes I doubt, since the output is usually 8 bits/color component.

nah i think 16 bits will be fine aswell. its going to be used to normalize 32bit/component normals. although the final output will be 8-bit colors, the intermediate output will be luminance values with a high range, but i dont think it will matter in the visual result anyway. i would like to experiment with it though, so i also want to take 16.16 roots in case in need them. its an interesting subject anyway, and i dont have any deadlines :).

Quote:

- Using SIMD technologies a sqrt would cost you a few portions of cycles. Not 200 !

a fixed point one? that sound really awesome.

Quote:

=> Thus, I'd not pay more than 1 cycle if it was for typical rendering.

what do you mean by that statement?

Quote:

Quote: Is this bad? I have the feeling it can be better. i dont want to resort to any float ops though. its quite time critical, because it will be used for normalization in a raytracer on a per-ray basis. ".. dont want to resort to any float ops though => its quite time critical ..." I hope that you did not mean an implication here. Because FPU is damn much speedier than integer arithmetics unless you consider SIMD instructions.

nono that was just a badly constructed sentence. im doing the fixed point thing mainly as a challenge, even if its twice as slow for some things, i dont care. those SIMD instructions sound interesting though, as long as it doesnt mean doing multiple roots in parralel. not that it isnt interesting, i just dont want to go over my head, i do like to finish this before i die :).

Quote:

There is far enough room for improvement if it's really "time critical". So be certain that your prerequisites are really valid and motivated.

good question. though it will get executed quite a few times, the speed of the application isnt likely to fall or stand even if its 200 cycles, although it will have an impact. im mainly just curious to know how im doing with my 65 cycles.

Quote:

Else you can also boost it using LUTs then NR. In fact it's how fast 3DNow or SSE rsqrt (thus sqrt=x*rsqrt) work internally. But the nature of floating point encoding (exponent and mantissa) makes the LUT approach very astute and efficient.

yeah i saw that method. i havnt benched it, but someone who seemed like he knew what he was talking about said it didnt work too good on modern processors because of all the branches. ill give it a try though.

Quote:

So once again, are you certain that you don't want any FPU or SIMD in your routine ? I don't mean necessarilly the sqrt of "math.h", the "noob" implementation.

nope, no fpu. id appriciate it if you could tell me a little more on SIMD and how it applies to squareroots though.

Quote:

Original post by Eelco nope, no fpu. id appriciate it if you could tell me a little more on SIMD and how it applies to squareroots though.

Most SIMD instruction sets have a hardware sqrt that is very fast if you're willing to sacrifice a little accuracy. For example, SSE has a 1 cycle reciprocal square root that is accurate to 12 bits.