/****************************************************
*   Program to demonstrate the Aitken Steffenson    *
*              iteration subroutine                 *
* ------------------------------------------------- *
* Reference: BASIC Scientific Subroutines, Vol. II  *
* By F.R. Ruckdeschel, BYTE/McGRAWW-HILL, 1981 [1]. *
*                                                   *
*               C++ version by J-P Moreau, Paris.   *
*                      (www.jpmoreau.fr)            *
* ------------------------------------------------- *
* SAMPLE RUN:                                       *
* (Example: find a real root of f(x) = x-2*SIN(x))  *
*                                                   *
* Input the initial guess:                          *
*                                                   *
*    X0 = -2                                        *
*                                                   *
* Convergence criterion: 1e-6                       *
*                                                   *
* Convergence factor: -1                            *
*                                                   *
* Maximum number of iterations: 25                  *
*                                                   *
*                                                   *
* The calculated zero is X = -1.895494              *
*                                                   *
* The associated Y value is Y = -0.000000           *
*                                                   *
* The number of iterations was: 3                   *
*                                                   *
* ------------------------------------------------- *
* Note: this algorithm fails with function (x+1)^5. *
****************************************************/
#include <stdio.h>
#include <math.h>

double  c,e,x,x0;
int     m, n;

//*************************************************
// Function subroutine                            
double Y(double x) {
  return (x-2.0*sin(x));
  c = 1.0-2.0*cos(x);
  c = -1.0/c;
}
//*************************************************


/**********************************************
*   Aitken Steffenson iteration subroutine    *
* ------------------------------------------- *
* This subroutine calculates the zeroes of a  *
* function Y(x) by iterations, and employs    *
* Aitken acceleration to speed up convergence.*
* An initial guess is required, x0, and two   *
* convergence factors, c and e. e relates to  *
* the accuracy of the estimate, and c is used *
* to aid convergence. Also required is the    *
* maximum number of iterations, m. c=-1 is a  *
* normal value, if divergence occurs, smaller *
* and/or positive values should be tried.     *
* The root is returned in x, the number of    *
* iterations in n.                            *
**********************************************/
void Steffenson()  {
  //Labels: e50,e100,e200
  double x1,x2,xk,yy;
  int m1;
  n=0;
e50: m1=0;
  x=x0;
  // Get yy
e100: yy = Y(x);
  yy=x+c*yy;
  // Enough points for acceleration ? 
  if (m1>0) goto e200;
  x1=yy; x=x1; n=n+1; m1=m1+1;
  goto e100;
e200: x2=yy;
  // Perform acceleration
  // Guard against a zero denominator
  xk=x2-2.0*x1+x0;
  if (xk==0) xk += 0.001;
  xk=(x1-x0)*(x1-x0)/xk;
  x0=x0-xk;
  // Test for convergence
  if (n>=m) return;
  if (fabs(x-x0)<e) return;
  x0=x1; x1=x2; x=x1;
  // Repeat process
  goto e50;
}


void main()  {
  printf("\n Input the initial guesses:\n\n");
  printf("    X0 = "); scanf("%lf",&x0);
  printf("\n Convergence criterion: "); scanf("%lf",&e);
  printf("\n Convergence factor: "); scanf("%lf",&c);
  printf("\n Maximum number of iterations: "); scanf("%d",&m);

  Steffenson();    // Call Aitken_Steffenson routine

  printf("\n\n The calculated zero is X = %f\n\n",x);
  printf(" The associated Y value is Y = %f\n\n",Y(x));
  printf(" The number of steps was: %d\n\n",n);
}

// End of file steffen.cpp