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

#define EPS 0.00000001
#define N   1000

typedef struct {
  double ans;
  int iter;
} Ans;

double f(double x)
{
  return 0.5 - x + 0.2 * sin(x);
}

/* f'(x) */
double f2(double x)
{
  return -1.0 + 0.2 * cos(x);
}

double g(double x)
{
  return 0.5 + 0.2 * sin(x);
}

Ans linear_iteration_method(double x, double (*f)(double))
{
  double x0;
  int i = 0;
  Ans ans;
  
  do {
    x0 = x;
    x = f(x);
    i++;
  } while (fabs(x - x0) > EPS);
  
  ans.ans = x;
  ans.iter = i;
  
  return ans;
}

Ans aitken_deltasquare_process(double a[], double b[])
{
  int i;
  double dx1, dx2;
  Ans ans;

  for (i = 0; i < N - 2; i++) {
    dx1 = a[i + 1] - a[i];
    dx2 = a[i + 2] - a[i + 1];
    b[i] = a[i + 2] - dx2 * dx2 / (dx2 - dx1);
    if (fabs(b[i] - a[i + 2]) <= EPS)
      break;
  }

  ans.ans = b[i];
  ans.iter = i + 1;

  return ans;
}

Ans newton_method(double x, double (*f)(double), double (*f2)(double))
{
  int i = 0;
  double x0;
  Ans ans;
  
  do {
    x0 = x;
    x = x - f(x) / f2(x);
    i++;
  } while (fabs(x0 - x) > EPS);
  
  ans.ans = x;
  ans.iter = i;
  
  return ans;
}

int main(void)
{
  int i;
  double x0 = 0.2;
  Ans ans;
  static double a[N], b[N];
  
  puts("線形反復法");
  ans = linear_iteration_method(x0, g);
  printf("solution = %.8f, iteration = %d\n", ans.ans, ans.iter);
  
  puts("エイトケンのデルタ2乗法");
  /* 元数列の生成 */
  a[0] = x0;
  for (i = 1; i < N; i++)
    a[i] = g(a[i - 1]);
  ans = aitken_deltasquare_process(a, b);
  printf("solution = %.8f, iteration = %d\n", ans.ans, ans.iter);

  puts("ニュートン法");
  ans = newton_method(x0, f, f2);
  printf("solution = %.8f, iteration = %d\n", ans.ans, ans.iter);

  return 0;
}
