0

Implement double sqrt(double x) in C++ without using std library.

This is a facebook interview question I saw here. http://www.glassdoor.com/Interview/Implement-double-sqrt-double-x-in-C-QTN_87210.htm Any other good idea about this?...

!!!Edited.!!!(without using std library.)

Josh Morrison
  • 7,488
  • 25
  • 67
  • 86

4 Answers4

8

Look here. This CodeProject article compares 14 different methods for computing the square root.

Lior Kogan
  • 19,919
  • 6
  • 53
  • 85
4

Here is one of the most genius sqrt implementations which can be found on wikipedia. It is not the most accurate, but extremely fast.

float fast_sqrt(float number) {
   long i;
   float x2, y;
   const float threehalfs = 1.5F;

   x2 = number * 0.5F;
   y  = number;
   i  = * ( long * ) &y;                     // floating point bit level hacking [sic]
   i  = 0x5f3759df - ( i >> 1 );             // Newton's approximation
   y  = * ( float * ) &i;
   y  = y * ( threehalfs - ( x2 * y * y ) ); // 1st iteration
   y  = y * ( threehalfs - ( x2 * y * y ) ); // 2nd iteration
   y  = y * ( threehalfs - ( x2 * y * y ) ); // 3rd iteration

   return 1/y;
}
Matthieu
  • 2,736
  • 4
  • 57
  • 87
regality
  • 6,496
  • 6
  • 29
  • 26
  • 6
    And you mean to tell us that's what you think the answer should have been at interview? – Tony Delroy Feb 15 '11 at 05:40
  • Did he want the job? You are supposed to try and impress them during the interview. – regality Feb 15 '11 at 06:18
  • @Tony: "I'd use the methodology of the fast inverse square root function from Quake - that formula is amazing - do you know if it?" would probably suffice. (though technically, this does 3 passes of newton-raphson approximation instead of 1... so it's not quite the same thing) – James Feb 15 '11 at 07:06
  • oh ya, I added those for a bit more accuracy, and forgot to take them out. basically the same thing tho. – regality Feb 15 '11 at 07:13
  • @regality: I've used this approach on integers and for 32bits a single iteration was sufficient for 100% accuracy (picked random numbers and including edges), have you checked the extra iterations for `float` ? – Matthieu M. Feb 15 '11 at 07:24
  • @James, @regality: spit it out from memory & it may well impress. But, what are they likely to expect? If it's not a job where you're that level of maths can be expected, or sqrt is particularly relevant, then just a working implementation of a practical, common-sense approach such as successive approximation. The fact that it's a constant power of 0.5 and not a log or other function is probably incidental. A reference to Quake would rate as a curiosity to me if you couldn't explain the algo, just as "C++ boost::any gives run-time typing" isn't worth much sans implementation insight. – Tony Delroy Feb 15 '11 at 07:31
  • OK I hope your okay with BoneOS using this.. We have cited its originality. Source https://github.com/Bone-Project/BoneOS/blob/amanuel_progress_kbd/libc/math/sqrtk/sqrtk.c – amanuel2 Nov 08 '16 at 16:57
4

The two obvious answers are bisection (semi-slow) and Newton-Raphson/Leibniz iteration (usually faster). To keep from spoiling anybody's fun, I'll do a reinterpret_cast on the question -- here's an implementation of an integer square root in 8086 assembly language using the Newton-Raphson technique:

isqrt proc uses di, number:word
;
; uses bx, cx, dx
;
    mov di,number
    mov ax,255
start_loop:
    mov bx,ax
    xor dx,dx
    mov ax,di
    div bx
    add ax,bx
    shr ax,1
    mov cx,ax
    sub cx,bx
    cmp cx,2
    ja  start_loop
    ret
isqrt endp

This is open to some improvement -- it uses x/2 as its initial guess at the sqrt(x). With 386+ instructions, you can use bsr to find the most significant bit that's set to get a rough approximation of log2x, and divide that by 2 to get your initial approximation.

OTOH, this really only made sense on ancient processors. For anything since the 486 (or so) that has built-in floating point hardware, it's nearly certain that the FSQRT instruction will beat this (or pretty much anything else you can write).

Jerry Coffin
  • 476,176
  • 80
  • 629
  • 1,111
  • as far as I know using the Reciprocal approach http://en.wikipedia.org/wiki/Methods_of_computing_square_roots#Reciprocal_of_the_square_root is faster, probably because it substitutes clever hacks to looping. I am unsure of its precision though. – Matthieu M. Feb 15 '11 at 07:26
  • @Matthieu M.: Yes, if you change the problem, you can improve speed a bit. Most of them also do a *very* crude approximation that only good for around 10 or 12 bits (which is still ~1 pixel at typical resolutions). For an actual square root, however, it's *really* hard to beat dedicated hardware. – Jerry Coffin Feb 15 '11 at 14:09
  • yes obviously, but I was referring to assembly function you wrote :) I've used the Reciprocal approach with perfect precision on integers, and I guess that adding some Newton's iterations would improve the precision somewhat. What's interesting is getting a good first approximation anyway, you can then refine the precision depending on your needs. – Matthieu M. Feb 15 '11 at 14:33
0

If I am allowed to use log (ln) and exp then of course exp(log(x)/2) will give me the square root.

Assuming not:

If our value we find the sqrt of is x and the start value is y then we iterate y->(y+x/y)/2

Terminating condition would either be a proximity of y to its previous value or of y*y to x.

With 385 as my x value I get these values in my iterations (Excel)

1
193
97.49740933
50.7231161
29.15667189
21.1805984
19.67880541
19.62150055
19.62141687
19.62141687

You can use an "approximate" 2^(log base 2(x)/2) as a start point instead of 1. 385 has a log somewhere between 8 and 9, so if we say 8.5 and therefore start with 2^4.25. If we do this linear between 16 and 32 then we would start with 20.

Starting with 20 I get there in just 4 steps:

20
19.625
19.6214172
19.62141687

but it required previous "iterations" to calculate the approximate log and exponential.

CashCow
  • 30,981
  • 5
  • 61
  • 92