#include <math.h>
#include <fstream.h>
#include <iomanip.h>
#include <vector>


// Chapter 4.1, p.149
// Polynomial interpolation
template <class Type>
void coef (int n, std::vector<Type> x, std::vector<Type> y, std::vector<Type>& a)
{
	int i, j;
	for (i=0; i<=n; i++)
		{ a[i] = y[i]; }
	for (j=1; j<=n; j++)
	{
		for (i=n; i>=j; i--)
			{ a[i] = (a[i] - a[i-1]) / (x[i] - x[i-j]); }
	}
}

template <class Type>
Type eval (int n, std::vector<Type> x, std::vector<Type> a, Type t)
{
	Type temp;
	temp = a[n];
	for (int i = n-1; i>=0; i--)
		{ temp = temp*(t-x[i]) + a[i]; }
	return (temp);
}


// Chapter 4.3, p.176-177
// Estimating Derivatives
template <class Type>
void deriv (Type f(Type), Type x, int n, Type h, std::vector< std::vector<Type> >& d)
{
	for (int i=0; i<=n; i++)
	{
		d[i][0] = (f(x+h) - f(x-h)) / (2.0*h);
		for (int j=0; j<=i-1; j++)
			{ d[i][j+1] = d[i][j] + (d[i][j] - d[i-1][j]) / (pow(4,j+1) - 1); }
		h *= 0.5;
	}
}