I(1)=1; J(50)=0; e=exp(1);
for n=2:50
    I(n)=e-n*I(n-1); % forward calculation
    k=50-n+2; J(k-1)=(e-J(k))/k;  % backward
end
fprintf('  n        I_forward            I_backward\n')
for n=1:50
    fprintf('%4i    %16.7e     %16.7e   \n',n,J(n),I(n))
end