%Lagrange Polynomial Algorithm by 
%Neville Interpolation at a vector of %points.




clc;
clear;
format long
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Input
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

n = input('Enter the number of points to interpolate at: ');
x = input('Enter the vector of points to interpolate at: ');
d = input('Enter the vector of function values of interpolation points: ');
t = input('Enter the value to approximate the function at: ');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Neville Method
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Q = zeros(n,n);	%Initializes an n x n matrix
Q(:,1) = d';	%Enters the vector of 
                % function values as the first column of Q.
Q;					

for i = 2:n
   for j = 2:i
      Q(i,j) =   ...
     [(t - x(i - j))*Q(i,j-1) - (t - x(i))*Q(i-1,j-1)]/(x(i) - x(i-j));
   end
end
%Q(i,j) calculates the polynomial approximation for each matrix point.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Output%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Q		%Displays the complete matrix Q.

P = input('Enter the point in the matrix to be outputed: ');
disp('The value of the polynomial is approximately: ')
disp(P)		%Displays the function approximation for that Lagrange
		%polynomial.