朗伯W函数

发布的克里夫硅藻土，2013年9月2日

77次(过去30天 0喜欢| 1评论

朗伯W函数应该更广为人知。它出现在各种各样的地方。我们用MATLAB函数求函数值是哈雷方法的一个很好的应用。

内容

约翰·兰伯特

约翰·亨利克·兰伯特于1728年出生在米卢斯，当时米卢斯在瑞士，但现在在法国。他是与瑞士著名数学家伦纳德·欧拉同时代的人。他的兴趣非常广泛。他是第一个证明$\pi$不合理的人。他引入了双曲三角函数。他对光学做出了如此重要的贡献，以至于光度单位以他的名字命名。他还在物理学、天文学、逻辑学和哲学方面做出了贡献。

兰伯特没有明确描述现在以他的名字命名的函数。他考虑了各种各样的方程，包括$x = q + x^m$这样的方程。欧拉接着又写了一些关于兰伯特工作的论文。分析这类方程所涉及的思想构成了我们今天所说的朗伯W函数。

枫连接

在20世纪90年代早期，Maple符号代数软件的开发者，包括Rob Corless、Gaston Gonnet、d.g. Hare和David Jeffrey，开始意识到这个函数在过去的几年里已经被重新发现了很多次。他们将这个函数命名为“W”，在枫树技术通讯上发表了一篇文章，并开始编写参考书目。Don Knuth通过贡献涉及树的枚举的应用程序加入了这个小组。1996年，他们五人写了一篇32页的论文，成为关于Lambert W函数的权威参考。他们的参考书目包含72条。

一个初等函数

让我们绘制一个初等函数，使用$w$作为自变量。函数是$w \mbox{e}^w$。

F = @(w) w.*exp(w);ezplot (f[4 1])包含(' w ')轴([-4 1 -。5 1])轴广场

你可以看到函数在$w = -1$处达到全局最小值。这个最小值在这个故事中扮演着重要的角色，所以让我们给它起个名字。

$ mbox{e} = -1/ mbox{e} $

Ebar = -1/exp(1)“标记”，“。”，“markersize”, 18岁,“颜色”，“黑”甘氨胆酸)组(,“ytick”, ebar 0、0.5、1)集(gca,“yticklabel”，{“1 / e”，' 0 '，“0.5”，' 1 '}）

ebar = -0.367879441171442

功能逆

朗伯W函数，$W(x)$，是$W \mbox{e}^ W $的逆函数。换句话说，$w$ = w ($x$)是方程的解

$ x = w \mbox{e}^w $

$W(x)$的图形在概念上是通过在上面的图形的垂直轴上加上an来获得的“x”然后交换水平轴和垂直轴。(你可以在下面的实际图中看到一个峰值，这是在我们看到如何计算函数后生成的。)

点$x = \bar{\mbox{e}} = -1/ mbox{e}$是至关重要的。在$\bar{\mbox{e}} <= x <= 0$区间上有两个实解，因此$W(x)$是双值的。金宝搏官方网站一个解是递增的，定义为所有$x >= \bar{\mbox{e}}$，是已知的主要分支，用$W_0(x)$表示。另一个解是递减的，它只定义为负$x$，极点为0，用$W_{-1}(x)$表示。(还有其他的复数分支，$W_k(x)$，这里我们没有考虑。)

哈雷的方法

哈雷法是一种求二阶导数为连续且容易计算的实值函数的零点的方法。它是以英国天文学家的名字命名的，因为他发现了一颗彗星。

牛顿的方法是用一个斜率相同的线性函数来近似一个函数，并且一步步逼近这个近似的零。哈雷的方法通过具有相同斜率和曲率的双曲线局部逼近函数，并逐步逼近该近似的最接近零的位置。迭代是

$ $ w_ {n + 1} = w_n-f左(w_n) / \ [f ' (w_n) - \压裂{(w_n) f”(w_n)} {2 f ' (w_n)} \] $ $

你可以通过忽略方括号中的第二项来理解牛顿方法。你可以看到，如果$f " (w)/f'(w)$恰好是一个简单的形式，哈雷的方法特别有效。但我们不认为哈雷法是一种通用方法，因为它需要二阶导数的知识。

朗伯W的应用

对于给定的$x$计算$W(x)$需要解这个方程

$ f(w) = w \mbox{e}^w - x = 0

w美元。我们可以用哈雷法，因为所需的导数是很容易得到的。

$ $ f ' (w) = w \ mbox {e} ^ w + w = \ \ mbox {e} ^ mbox {e} ^ w (w + 1) $ $

$ $ f”(w) = w \ mbox {e} ^ w + w \ mbox {e} ^ w + w = \ \ mbox {e} ^ mbox {e} ^ w (w + 2) $ $

请注意,

$ f " (w)/f'(w) = (w+2)/(w+1) $

这给了我们在给定值$x$下计算Lambert W的以下迭代步骤。

$ f_n = w_n \mbox{e}^{w_n} - x $$

$ $ w_ {n + 1} = w_n-f_n /左\ [\ mbox {e} ^ {w_n} (w_n + 1) - (w_n + 2) fn / (2 w_n + 2) \] $ $

起始值

在我之前，已经有很多人编写了程序，用哈雷彗星的方法来计算Lambert W。他们对减少最终迭代次数的初始值进行了相当多的关注。我决定采取一种不同的方法。如果程序能简化，我愿意把计算再进行几次迭代。事实上，我用的是最简单的初值。

对于主分支，$x$从$w = 1$开始。不要试图接近$W(x)$。对于任何$x$，迭代结果在大约6次迭代中收敛到双精度浮点精度，除非值非常接近分支点$\bar{\mbox{e}}$。更准确的起始值可以将迭代次数减半，但我更喜欢更简单的代码。

对于较低的分支，从任何小于$-1$的值开始。我选择了$w = -2$，但我不认为这很重要。同样，除非$x$非常接近分支点$\bar{\mbox{e}}$或极点$0$，否则大约需要六次迭代。

停止值

哈雷法在单根附近是三次收敛的，就是这个。一旦我们接近了，我们就会迅速接近。如果两个连续迭代一致1.0 e-8，那么你可以相当肯定第二个是准确的1.0 e-16．我甚至可以换一个1. e-8通过eps ^ (1/3) = 6 e-6．

Lambert_W.m

这是整个m文件函数Lambert_W．这将在适当的时候从MATLAB中央文件交换可用。(添加9/16/2013:明白了这个链接）

类型Lambert_W

函数的逆(x = w*exp(w))。% w = Lambert_W (x)一样Lambert_W (x, 0) % w = Lambert_W (x, 0)主键或上分支,W_0 (x) % w = Lambert_W (x, 1)较低的分支,W_ {1} (x) % %见:https://blogs.mathworks.com/cleve/2013/09/02/the-lambert-w-function/ % 2013年版权MathWorks公司%有效开始猜如果输入参数个数< 2 | | ~ = 1 w =分支的(大小(x));% Start above -1 else w = -2*ones(size(x));% Start below -1 end v = inf*w;% Haley's method while any(abs(w - v)./abs(w) > 1.e-8) v = w;e = exp (w);F = w.*e - x;%迭代使这个量0 w = w - f / ((e。* (w + 1)——(w + 2)。* f / w + (2 * 2)));结束

图兰伯特W

我们现在可以用这段代码来绘制朗伯W函数的两个分支的图。

%的主要分支x = ebar: .01:1;情节(x, Lambert_W (x, 0))%降低分支持有在x = ebar: .01:0;情节(x, Lambert_W (x, 1),“r”)举行从%注释轴([-。51-4 1]) axis广场线([ebar ebar], [1],“标记”，“。”，“markersize”, 18岁,“颜色”，“黑”甘氨胆酸)组(,“xtick”，[ebar 0 .5 1])“xticklabel”，{“1 / e”，' 0 '，“0.5”，' 1 '})集(gca,“ytick”4:1)包含(“x”) ylabel (' w ')标题(“兰伯特W”)({传奇“W_0”，“W_{1}”}，“位置”，“东南”）