sqrt.c (675B)
1 /* 2 sqrt returns the square root of its floating 3 point argument. Newton's method. 4 5 calls frexp 6 */ 7 8 #include <u.h> 9 #include <libc.h> 10 11 double 12 sqrt(double arg) 13 { 14 double x, temp; 15 int exp, i; 16 17 if(arg <= 0) { 18 if(arg < 0) 19 return 0.; 20 return 0; 21 } 22 x = frexp(arg, &exp); 23 while(x < 0.5) { 24 x *= 2; 25 exp--; 26 } 27 /* 28 * NOTE 29 * this wont work on 1's comp 30 */ 31 if(exp & 1) { 32 x *= 2; 33 exp--; 34 } 35 temp = 0.5 * (1.0+x); 36 37 while(exp > 60) { 38 temp *= (1L<<30); 39 exp -= 60; 40 } 41 while(exp < -60) { 42 temp /= (1L<<30); 43 exp += 60; 44 } 45 if(exp >= 0) 46 temp *= 1L << (exp/2); 47 else 48 temp /= 1L << (-exp/2); 49 for(i=0; i<=4; i++) 50 temp = 0.5*(temp + arg/temp); 51 return temp; 52 }