
function b = wynneps(a);

%  Quick and dirty implementation of Wynn's epsilon algorithm.
%  The number b is the computed limit of the sequence a.
%  Warning:  Nonlinear sequence accelerators are known to be unstable.
%  No special provision to avoid this was made in the code below; 
%  use at your own risk.

%                                 Andre Weideman, Updated: November 1997

a = a(:);                               % Make sure a is a column vector.
N = length(a);

E(:,1:2) = [zeros(size(a)) a];

for j = 3:N
d = [E(2:N,j-1); 0]-[E(1:N-1,j-1); 0];  % This finagling avoids Matlab's
J = (d==0);  d = d+eps*J;               % division by zero error messages.  
E(:,j) = [E(2:N,j-2); 0] + 1./d;        
end

b = E(1,2:2:N); b = b(length(b));

















