pracma (version 1.7.0)

lambertWp: Lambert's W Function

Description

Principal branch of the Lambert W function

Usage

lambertWp(z)
D_lambertWp(z)

Arguments

z

Numeric vector of real numbers >= -1/e.

Value

• Returns the solution w of w*exp(w) = z for real z with NA if z < 1/exp(1).

Details

The Lambert W function is the inverse of x --> x e^x, which is unique for x >= -1/e. Here only the principal branch is computed for real z.

The value is calculated using an iteration that stems from applying Newton's method. This iteration is quite fast.

D_lambertWp returns the exact derivative of lambertWp, based on the fact that lambertWp is the inverse of x*exp(x).

The functions is not really vectorized, but at least return a vector of values when presented with a numeric vector of length >= 2.

References

Corless, R. M., G. H.Gonnet, D. E. G Hare, D. J. Jeffrey, and D. E. Knuth (1996). On the Lambert W Function. Advances in Computational Mathematics, Vol. 5, pp. 329-359. http://www.apmaths.uwo.ca/~djeffrey/Offprints/W-adv-cm.pdf.

See Also

newtonRaphson

Examples

##  Examples
lambertWp(0)          #=> 0
lambertWp(1)          #=> 0.5671432904097838...  Omega constant
lambertWp(exp(1))     #=> 1
lambertWp(-log(2)/2)  #=> -log(2)

# The solution of  x * a^x = z  is  W(log(a)*z)/log(a)
# x * 123^(x-1) = 3
lambertWp(3*123*log(123))/log(123)  #=> 1.19183018...

xs <- c(-1/exp(1), seq(-0.35, 6, by=0.05))
ys <- lambertWp(xs)
plot(xs, ys, type="l", col="darkred", lwd=2, ylim=c(-2,2),
     main="Lambert W0 Function", xlab="", ylab="")
grid()
points(c(-1/exp(1), 0, 1, exp(1)), c(-1, 0, lambertWp(1), 1))
text(1.8, 0.5, "Omega constant")

# Second branch resp. the complex function lambertWm()
F <- function(xy, z0) {
    z <- xy[1] + xy[2]*1i
    fz <- z * exp(z) - z0
    return(c(Re(fz), Im(fz)))
}
newtonsys(F, c(-1, -1), z0 = -0.1)   #=> -3.5771520639573
newtonsys(F, c(-1, -1), z0 = -pi/2)  #=> -1.5707963267949i = -pi/2 * 1i