# Compute the Lambert W function

Your challenge is to compute the Lambert W function. The W of x is defined to be the real value(s) y such that

y = W(x) if x = y*(e^y)


where e = 2.718281828 is Euler's number.

Sometimes, y may not be real.

### Examples

W(-1) = non-real
W(-0.1) = -0.11183, -3.57715
W(1) = 0.56714
W(2) = 0.85261


Here's a quick graph of what this function looks like.

### Rules

Your goal is to take an input and output either nothing, 1 solution, or 2 solutions, out to 5 significant figs. You should expect float inputs within the reasonable range of -100..100.

This is code-golf, so shortest code wins.

• A note about Lambert W function: there is no solution for x < -1/e, two solutions for -1/e < x < 0, and one solution for x == -1/e or x >= 0. – JungHwan Min Nov 23 '17 at 5:58
• 5 places after the decimal, or 5 significant figures? Since W^-1(x) tends to negative infinity when x->0-, it is hard to get 5 places after the decimal. – Colera Su Nov 23 '17 at 8:04
• If one input should give two values, this is not a function – Luis Mendo Nov 23 '17 at 10:43
• @Mr.Xcoder Sure, though the non-real output should be reasonable. – Simply Beautiful Art Nov 23 '17 at 12:50
• I think a nice way to handle multiple solutions is the output being a list of reals, containing 0, 1, or 2 values. Would this be allowed? – xnor Nov 23 '17 at 16:56

# Octave, 36 34 bytes

@(x)(y=lambertw(-1:0,x))(~imag(y))


Try it online!

I solved this in 10 different ways without the builtin, but I couldn't get the solution for both the -1 branch and the 0-branch. I wasn't able to solve this without resorting to the builtin lambertw function. Even though it uses a builtin, it's far from straightforward.

lambertw will only give out one value by default, namely the W0(x), unless otherwise specified. The problem arises when we specifies that we want two values, both W0(x) and W-1(x). If this is the case, and x>=0, it will return one real and one complex result.

We therefore have to specify that we only want the results where the imaginary part of the result is zero.

### Breakdown:

@(x)                              % Anonymous function that takes x as input
(y=lambertw(-1:0,x))          % The lambertw function, with both the -1 and 0 branches
(~imag(y))  % Index, where the elements with imag(y)=0 is true


# 05AB1E, 98867574 73 bytes

-11 bytes thanks to @KevinCruijssen

žr(zV"Ðžrsm*I-d"UY‹iõë∞.ΔX.V}®s‚ΔDÅAX.Vǝ}н,0‹i∞.Δ(X.V}(®‚ΔDÅAX.V≠ǝ}IYÊiн,


Try it online!

Explanation:

(I will let f(x)=x*e^x)

žr(zV"Ðžrsm*I-d"U #Let Y=-1/e,X=str(ÐYsm*I-d) to save characters in the future

Y‹iõ              #if input<-1/e: print empty string

ë∞.ΔX.V}          #else: find least int k>=1 such that f(k)>=0 (evaluates string X as though it were a function)

®s‚               #make list [-1,k]

ΔDÅAX.Vǝ}         #find root of f(x)-input in the interval [-1,k] using bisection method

#actually this outputs a list with 2 numbers which will be very close to each other, which represents the interval in which the root lies

#e.g. the list [-3.00001,-3] would represent the interval [-3.00001,-3]

н,                #print the first item of the list; stack is currently empty

0‹i               #if input<0:

∞.Δ(X.V}(         #find greatest int k<=-1 such that f(k)>=0

®‚                #make list [k,-1]

ΔDÅAX.V≠ǝ}        #find root of f(x)-input in the interval [k,-1] using bisection method

#again, this outputs a list rather than a single number

IYÊiн,            #if input!=-1/e, print that root (done so that "-1" is not printed twice)


Visualization of the program's indentation:

žr(zV"Ðžrsm*I-d"UY‹iõë∞.ΔX.V}®s‚ΔDÅAX.Vǝ}н,0‹i∞.Δ(X.V}(®‚ΔDÅAX.V≠ǝ}IYÊiн,
{    {   }   {       }    {  {    }   {        }   {

• I'm sure this can be golfed some more, but here is 75 bytes by replacing both 1( with ® (relevant 05AB1E tip), and putting ÐYsm*I-d inside a string variable at the start, so we can use it four times with .V. – Kevin Cruijssen Apr 20 at 12:10

# Ruby, 131 bytes

->x{y,z,e,l=1,-1,Math::E,->t{Math::log t};100.times{z=x<0&&x>-1/e ?l[x/z]:p;y=x<-1/e ?"non-real":x>e ?l[x/y]:x*e**-y};p y;z&&(p z)}


Try it online!

Basically applies fixed-point iteration 100 times, as partially explained in Lambert W function aproximation.

• Why non-competing? – MD XF Nov 23 '17 at 3:38
• Marking a submission non-competing does not make it valid. If your submission is invalid, please edit it to make it valid or delete it altogether. Relevant meta post. – JungHwan Min Nov 23 '17 at 4:51
• Actually the only non-competing reason I have seen is language postdate the challenge. In other cases I often see they are deleted or competing. / Well, apply the fixed-point iteration is a valid way to compute the function, but only if you can prove that it is good enough (to your own requirement). – user202729 Nov 23 '17 at 9:31
• @JungHwanMin My submission is valid, I just marked it non-competing because, well, I didn't want to compete. I suppose that in itself is not a valid reason :P – Simply Beautiful Art Nov 23 '17 at 12:49

# C, 101 + 2 (-lm) bytes

Returns how many real solutions, and store them in an array.

This one is very slow, since it uses a brute-force approach to find roots. It handles up to x = 9e9.

f(a,b,p,t,c,x)float*b,a,x;{for(p=a<0,c=0,x=-50;x<9;x+=2e-6,p=t)(t=x*exp(x)>a)^p?b[c++]=x:0;return c;}


Try it online!

# C, 148 146 + 2 (-lm) bytes

This function returns how many real solutions, and store them in the given pointers.

Use binary search to find a solution. It handles up to about x = 108.

#define h(x)99;for(l=-1;x=(l+r)/2,fabs(l-r)>1e-6;)x*exp(x)>a?r=x:(l=x);
f(a,b,c,l,r)float*b,*c,a,l,r;{r=h(*b)r=-h(*c)return-1/exp(1)>a?0:a<0?2:1;}


Try it online!

# Wolfram Language (Mathematica), 28 bytes

x/.NSolve[#==E^x x,x,Reals]&


Try it online!

Returns x for no-solution.

This function should give all solutions because NSolve algorithm uses inverse functions, and there is a built-in inverse function of x = y*(e^y) called ProductLog.

### Explanation

#==E^x x


Equation <input> = e^x * x

NSolve[ ... ,x,Reals]


Find solutions in the domain of real numbers.

x/.


Replace x with the solutions.

### Alternative version, 35 bytes

N@ProductLog[{-1,0},#]~Cases~_Real&


Try it online!

# Python 3, 258 bytes

e,f,X,R,L=2.7182818,lambda x:x*e**x,float(input()),0,-1
if X>=-1/e:
while f(R)<X:R+=1
while f(R)>X:
M=(L+R)/2
if f(M)<X:L=M
else:R=M
print(M)
if X<0:
while f(L)<X:L-=1
while f(L)!=X:
M=(L+R)/2
if f(M)<X:R=M
else:L=M
if M!=-1:print(M)


Try it online!

A port of my 05AB1E answer

# Casio BASIC, 192 Bytes

A port of My 05AB1E answer

?→X
""
If X≥-e^-1
Then 0→R
While Re^R<X
Isz R
WhileEnd
-1→L
While Re^R>X
.5(L+R→M
If Me^M<X
Then M→L
Else M→R
IfEnd
WhileEnd
If X<0
Then M◢
While Le^L<X
Dsz L
WhileEnd
While Le^L≠X
.5(L+R→M
If Me^M<X
Then M→R
Else M→L
IfEnd
WhileEnd
If M≠-1
Then M


Note that indents are actually not supported; I added them for legibility.