#include <math.h>
#include <stdio.h>
#include <time.h>


#define FLOOR 1.0f
#define CEILING 100000.0f
#define STEP 100

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

	x2 = number * 0.5F;
	y  = number;
	i  = * ( long * ) &y;						
	i  = 0x5f3759df - ( i >> 1 );               
	y  = * ( float * ) &i;
	y  = y * ( threehalfs - ( x2 * y * y ) );   // 1st iteration
//	y  = y * ( threehalfs - ( x2 * y * y ) );   // 2nd iteration, this can be removed

	return y;
}


float simple_rsqrt(float x) {
    return 1.0f / sqrtf(x);
}



double simple_rsqrt_benchmark() {
    float x, y;
    x = FLOOR;
    clock_t begin = clock();
    while (x < CEILING) {
        x = nextafterf(x, CEILING);
        y = simple_rsqrt(x);
    }
    clock_t end = clock();
    return (double)(end - begin) / CLOCKS_PER_SEC;
 
}


double Q_rsqrt_benchmark() {
    float x, y;
    x = FLOOR;
    clock_t begin = clock();
    while (x < CEILING) {
        x = nextafterf(x, CEILING);
        y = Q_rsqrt(x);
    }
    clock_t end = clock();
    return (double)(end - begin) / CLOCKS_PER_SEC;
 
}

int bench_diff() {
    float x, y;
    x = FLOOR;
    while (x < CEILING) {
        for (int i = 0; i < STEP; i++) {
            x = nextafterf(x, CEILING);
        }
        y = fabsf(Q_rsqrt(x) - simple_rsqrt(x));
        printf("%.8f %.8f\n", x, y);
    }
    return 0; 
}




int main() {
    double x = simple_rsqrt_benchmark();
    double y = Q_rsqrt_benchmark();
    printf("Time for simple rqsrt %g\n", x);
    printf("Time for quakes rqsrt %g\n", y);
    return 0;
}
