clc,clear

%% 初始化
A = zeros(3,4);
A(:,1) = [0.5,0.6,0.4];
A(:,2) = [-0.6931,-0.5108,-0.9163];
syms x

%% 计算差商
for m = 3:size(A,2)
    
    A(m-1:end,m) = diff( A(m-2:end,m-1) )./diff( A(1:m-2:end,1) ); 
    
end

%% 牛顿插值多项式: for循环
% N = A(1,2);
% omega =1; 
% for n = 2:size(A,1)
%     
%     omega = omega*(x-A(n-1,1));
%         N = N + A(n,n+1)*omega;
%     
% end

%% 牛顿插值多项式: 无for循环
B = [ 1; x - A(1:end-1,1) ];
C = B * ones( size(B') );
D = triu(ones(size(C)),0);
E = D.*C + tril(ones(size(C)),-1);
N = prod(E) * diag(A(:,2:end));


%% 输出
fprintf('\t x \tf(x) \t一阶差商   二阶差商\n')
disp(A);
fprintf('N(x)=%s\n',vpa(simplify(N)));
fprintf('N(0.54)=%s\n',vpa(subs(N,0.54)));