(28/48/50) Lambert W Function

[Albert Chan]

(03-31-2023 10:07 PM)Albert Chan Wrote:  We may not need accurate log1p(x)-x after all.

With accurate log1p(x)-x, we can get W around branch point to within 1 ULP or less.
With below FNL(X), previous post examples all match W exactly, no ULP errors.

200 DEF FNL(X) ! = ln(1+X) - X, but more accurate
210 IF ABS(X)>=.4 THEN X=X/(SQRT(1+X)+1) @ FNL=FNL(X)*2-X*X @ GOTO 250
220 X2=X/(X+2) @ X4=X2*X2
230 X4=X4*(5005-X4*(5082-X4*969))/(15015-X4*(24255-X4*(11025-X4*1225)))
240 FNL=X2*(X4+X4-X)
250 END DEF

line 230: X4 = pade approx. of atanh(X2)/X2 - 1 ≈ ln(1+X)/(2*X2) - 1
line 240: FNL = X2*(X4+X4-X) ≈ ln(1+X) - 2*X2 - X*X2 = ln(1+X) - X

42 REPEAT @ L=FNL(X) @ Z=X-(H-L)/(X+L) @ X=X-Z @ UNTIL X=X+Z*.0001
44 PRINT "e^W, W =";R+R*X;LOGP1(X)-1

Away from branch point is a different story. With X getting big, FNL(X) won't help much.