/*         ATK IV - Numerical Programming (2006)         */
/* Excercise 4.4 - Rational and Polynomial Interpolation */

#include<stdio.h>
#include<math.h>
#define NRANSI
#include<nr.h>

#define N 20
#define Pi 3.1415962

float f(float x){
  return log(x*x/2.0);
}

int main(void){
  float xa[N+1],ya[N+1];
  float pol_y, rat_y, dy;
  float x;
  float xi = 0.0000001;
  float xf = 1.0/10.0;
  float step = (xf-xi)/(1.0*(N-1));
  int i;
  
  x=xi;
  for(i=1;i<=N;i++){
    xa[i] = x;
    ya[i] = f(x);
    /* printf("%f %f \n",x,f(x));*/
    x += step;
  }

  x=xi;
  step = (xf-xi)/(1.0*(10*N-1));
  for(i=0;i<10*N;i++){
    polint(xa,ya,N,x,&pol_y,&dy);
    ratint(xa,ya,N,x,&rat_y,&dy);
    printf("%f %f %f %f\n",x,pol_y,rat_y,f(x));
    x+=step;
  }
  return 0;
}