x_nodes = [1 4 9 16 25];
y_values = [1 2 3 4 5];
t = 0:.01:max(x_nodes);



for i = 2:length(x_nodes)
    
divdif_y = divdif(x_nodes(1:i),y_values(1:i));

p_eval = interp(x_nodes(1:i),divdif_y,t);


figure(i-1)
axis([1 25 1 6])
hold all
plot(t,0)
plot(t,p_eval,'Linewidth',2,'Color','red')
plot(x_nodes,y_values,'*','Color','black')
end