/*  Calcolo iterativo di radice quadrata di numeri positivi

    1)  se x > sqrt(N) ->        x^2 > N 
                           x^2 + x^2 > N + x^2
                              2*x^2  > N + x^2
                                 2*x > N/x + x 
                                   x > (x + N/x) / 2

      Ovvero, se x > sqrt(N), x è maggiore della media aritmetica 
                              fra x stesso e N/x 

    2) Ma una media aritmetica di due numeri positivi è sempre 
       maggiore della media geometrica, ovvero
                   (x + N/x) / 2  > sqrt( x * N/x ) = sqrt(N) 


    1)+2):   sqrt(N) < (x + N/x) / 2 < x
                               
    3) Ma, se chiamiamo x1 = (x + N/x) / 2, essendo > sqrt(N),
       segue

             sqrt(N) < (x1 + N/x1) / 2 < x1

    Etc. etc.  

      => si può cominciare con x=N

 */


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

int main(int argc, char *argv[]) {
  double x, N;
  int i;

  if (argc < 2) {
    printf("Sintassi:  ./sqrt_babilonesi N\n");
    return 0;
  }
  N = atof(argv[1]);

  x = N;                             // ma si può anche
                                     // partire da un numero >> N  (!!)
  printf("i=%2d: %.17f\n", 0, x);
  for (i=1; i<=10; i++) {
    x = (x+N/x) / 2.;
    printf("i=%2d: %.17f\n", i, x);    
  }
  printf("\nCon sqrt(): %.17f\n", sqrt(N));

  printf("\nCon bc:     1.41421356237309504880\n");
  printf("\[ echo 'sqrt(2)' |  bc -l ]\n");
  
  return 0;
}