samer@4: function W = lambertw(x)
samer@4: % lambertw - Lambert's W function
samer@4: %
samer@4: % lambertw :: real -> real.
samer@4: %
samer@4: % evaluates W(y) for -1 branch of Lambert W function where
samer@4: % y = -exp(-x)
samer@4: %
samer@4: % The W function satisfies W(t)*exp(W(t)) = t
samer@4: 
samer@4: if x<1, error('lambertw: argument must be >= 1'); end
samer@4: W=fixpoint(@(w)(-x-log(abs(w))),-x,'its',60);