#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>

double sqrt_knoppi(double a) {
    double xneu = 0, xalt = 0, epsilon = 0.00000000001;
 
    xneu = a;
    do {
        xalt = xneu;
        xneu = 0.5*(xalt+a/xalt);
    } while((xalt-xneu) > epsilon);
 
    return xneu;
}

const int64_t magic_number = 0x5fe6eb50c7b537a9;
double  sqrt_q3(const double x)
{
  const double xhalf = 0.5f*x; 
  union
  {
    double x;
    int64_t i;
  } u;
  u.x = x;
  u.i = magic_number - (u.i >> 1); 
  return x*u.x*(1.5f - xhalf*u.x*u.x);
}   

double  sqrt_q3_2(const double x)
{
  const double xhalf = 0.5f*x; 
  union
  {
    double x;
    int64_t i;
  } u;
  u.x = x;
  u.i = magic_number - (u.i >> 1); 
  u.x = u.x*(1.5f - xhalf*u.x*u.x); // Mehr von dieser Zeile für höhere Präzision
  return x*u.x*(1.5f - xhalf*u.x*u.x);
}   

void print_diff(double result, double sqrt)
{
  printf("%.15f \tdiff = %.15f\tRel = %.15f\n", result, fabs(sqrt -result), fabs(sqrt -result)/sqrt);
}

void sqrt_test(double d)
{
  double result = sqrt(d);
  printf("sqrt(%.15f) = %.15f\n", d, result);
  printf("knoppi: sqrt(%.15f) = ", d); print_diff(sqrt_knoppi(d), result);
  printf("knoppi: sqrt_q3(%.15f) = ", d); print_diff(sqrt_q3(d), result);
  printf("knoppi: sqrt_q3_2(%.15f) = ", d); print_diff(sqrt_q3_2(d), result);
  putchar('\n');
}

int main()
{
  double d = 15;
  for(int i = 0; i < 30; ++i)
    {
      sqrt_test(d);
      d /= 3;      
    }
}
