// rho computation with floats for homework, by Frank Ruskey. // Method 1 uses "uint64_t" from stdint; could also use "unsigned long long". // Method 2 uses a union from idea of Riley Mottus. #include #include typedef union { double d; struct { uint64_t mant : 52; uint16_t expo : 11; uint16_t sign : 1; } parts; } dull; int main( int argc, char *argv[] ){ uint64_t x, p, two20=1ULL<<20; for (x=two20+1; x<=two20+32; ++x) { // Method 1 double f64 = (double)(x & -x); p = ( *(( (uint64_t *) &f64 ))>> 52 ) - 1023ULL; // Method 2 dull y; y.d = (double)(x & -x); printf( "rho(%llu) = %llu %d\n", x, p, y.parts.expo-1023 ); } return(0); } /* Frank@Frank-PC ~ $ gcc RhoTest.c Frank@Frank-PC ~ $ ./a.exe rho(1048577) = 0 0 rho(1048578) = 1 1 rho(1048579) = 0 0 rho(1048580) = 2 2 rho(1048581) = 0 0 rho(1048582) = 1 1 rho(1048583) = 0 0 rho(1048584) = 3 3 rho(1048585) = 0 0 rho(1048586) = 1 1 rho(1048587) = 0 0 rho(1048588) = 2 2 rho(1048589) = 0 0 rho(1048590) = 1 1 rho(1048591) = 0 0 rho(1048592) = 4 4 rho(1048593) = 0 0 rho(1048594) = 1 1 rho(1048595) = 0 0 rho(1048596) = 2 2 rho(1048597) = 0 0 rho(1048598) = 1 1 rho(1048599) = 0 0 rho(1048600) = 3 3 rho(1048601) = 0 0 rho(1048602) = 1 1 rho(1048603) = 0 0 rho(1048604) = 2 2 rho(1048605) = 0 0 rho(1048606) = 1 1 rho(1048607) = 0 0 rho(1048608) = 5 5 */