/* [wxMaxima batch file version 1] [ DO NOT EDIT BY HAND! ]*/
/* [ Created with wxMaxima version 0.8.1 ] */

/* [wxMaxima: title   start ]
Pi - according to the method of Archimedes 
(in the form of Christian Wolff) 
   [wxMaxima: title   end   ] */

/* [wxMaxima: comment start ]
Archimedes of Syracuse, approx. 287-212 BC 
Christian Freiherr von Wolff, approx. 1679-1754 AD 
   [wxMaxima: comment end   ] */

/* [wxMaxima: comment start ]
Reference: J. Ziegenbalg: Algorithmen von Hammurapi bis Gödel,  
  Verlag Harri Deutsch, Frankfurt a.M. 2007, section 3.2.4 
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
Pi_Archimedes_Wolff(s) := 
  block([r:1, se, su, ue, uu, i, n:3], 
    se : sqrt(3),  /* Belegung mit den Anfangswerten */    
    ue : 3 * se,   /* d.h.:  Werte fuer das Dreieck  */
    su : 2 * sqrt(3), 
    uu : 3 * su, 
    print(n, "  ", ue/2, "  ", uu/2, "  ", se*se), 

    for i : 1 step 1 thru s do 
      (n  : n * 2, 
       se : r*sqrt(2-2*sqrt(1-(se/(2*r))*(se/(2*r)))), 
       ue : n * se, 
       su : se / sqrt(1 - (se/(2*r)) * (se/(2*r)) ), 
       uu : n * su, 
       print(i, n, ue/2, uu/2, se*se) ),  
    [(ue/2+uu/2)/2, bfloat((ue/2+uu/2)/2)] );
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
Pi_Archimedes_Wolff(5);
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
fpprec : 100 $
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
set_display(ascii);
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
Pi_Archimedes_Wolff_printing_as_bigfloat(s) := 
  /* in this version all printing is done in the bigfloat mode */ 
  block([r:1, se, su, ue, uu, i, n:3], 
    se : sqrt(3),  /* initial values             */    
    ue : 3 * se,   /* for the "triangle"-polygon */
    su : 2 * sqrt(3), 
    uu : 3 * su, 
    printf(true, "~2d ~10d ~13,10h ~13,10h ~43,40h ~%", 
                      0, n, bfloat(ue/2), uu/2, se*se), 
    for i : 1 step 1 thru s do 
      (n  : n * 2, 
       se : r*sqrt(2-2*sqrt(1-(se/(2*r))*(se/(2*r)))), 
       ue : n * se, 
       su : se / sqrt(1 - (se/(2*r)) * (se/(2*r)) ), 
       uu : n * su, 
       printf(true, "~2d ~10d ~13,10h ~13,10h ~43,40h ~%", 
                      i, n, ue/2, uu/2, se*se) ), 
    bfloat((ue/2+uu/2)/2) );
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
Pi_Archimedes_Wolff_printing_as_bigfloat(30);
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
 Pi_Archimedes_Wolff_computation_as_bigfloat(s) := 
  /* in this version all computation is done in the bigfloat mode */ 
  block([r:1, se, su, ue, uu, i, n:3], 
    se : bfloat(sqrt(3)),  /* initial values             */    
    ue : bfloat(3 * se),   /* for the "triangle"-polygon */
    su : bfloat(2 * sqrt(3)), 
    uu : bfloat(3 * su), 
    printf(true, "~2d ~10d ~13,10h ~13,10h ~43,40h ~%", 
                      0, n, ue/2, uu/2, se*se),
    for i : 1 step 1 thru s do 
      (n  : bfloat(n * 2), 
       se : bfloat(r*sqrt(2-2*sqrt(1-(se/(2*r))*(se/(2*r))))), 
       ue : bfloat(n * se), 
       su : bfloat(se / sqrt(1 - (se/(2*r)) * (se/(2*r)))), 
       uu : bfloat(n * su), 
       printf(true, "~2d ~10d ~13,10h ~13,10h ~43,40h ~%", 
                      i, n, ue/2, uu/2, se*se) ), 
    bfloat(0.5*(0.5*ue+0.5*uu)) );
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
Pi_Archimedes_Wolff_computation_as_bigfloat(30);
/* [wxMaxima: input   end   ] */

/* Maxima can't load/batch files which end with a comment! */
"Created with wxMaxima"$