27 lines
532 B
C
27 lines
532 B
C
|
#include <math.h>
|
||
|
|
|
||
|
|
typedef union {
|
||
|
|
double d;
|
||
|
|
struct {
|
||
|
|
unsigned long int mantissa : 52;
|
||
|
|
unsigned int exponent : 11;
|
||
|
|
unsigned int sign_bit : 1;
|
||
|
|
} parts;
|
||
|
|
} _dbl_cast;
|
||
|
|
|
||
|
|
double sqrt(double x) {
|
||
|
|
_dbl_cast S = { .d = x };
|
||
|
|
double estimate = (0.5 + 0.5 * S.parts.mantissa) *
|
||
|
|
exp2(S.parts.exponent - 1023);
|
||
|
|
|
||
|
|
double a_n, b_n;
|
||
|
|
|
||
|
|
for(int i = 0; i < 15; i++) {
|
||
|
|
a_n = (x - (estimate * estimate)) / (2 * estimate);
|
||
|
|
b_n = estimate + a_n;
|
||
|
|
estimate = b_n - ((a_n * a_n) / (2 * b_n));
|
||
|
|
}
|
||
|
|
|
||
|
|
return estimate;
|
||
|
|
}
|