math - 求解包含二 Gamma 函数的方程组的最有效方法是什么？

Question

求解涉及 digamma 函数的方程组的最有效方法是什么？

我有一个向量 v，我想求解一个向量 w，使得对于所有 i:

digamma(sum(w)) - digamma(w_i) = v_i

和

w_i > 0

我找到了 gsl 函数 gsl_sf_psi，它是 digamma 函数(使用某种级数计算)。是否有一个恒等式可以用来简化方程？我最好的选择是使用求解器吗？我正在使用C++0x；哪种求解器最容易使用且速度最快？

<小时/>

根据我的初步研究，digamma 不容易可逆(搜索逆 digamma 给出了通过二分搜索工作的算法)，因此整个系统没有简化是有道理的。

因此，使用求解器现在会留下两个问题:处理 digamma 计算速度非常慢的事实，以及处理 w_i > 0 的限制，否则 digamma(w_i) 将在 w_i = 0 时崩溃。

对于第一个问题，我想也许我应该为最近计算的 digamma 值实现一个缓存 - 我认为这是一个好主意，但不太了解寻根算法的工作原理。

我的想法是解决第二个问题是找到w'_i = log(w_i)。这样，w'_i就在整条线上。我想知道这是否是一个好主意。可能没有直接找到 digamma(exp(w')) 的函数？此外，该算法可能在 w' 空间中采取步骤，但不会改进事情，因为从 w'->w 的映射会失去一些精度，因此 w' 的两个元素可能映射到相同的 w。

仍然存在找到一个好的、快速的寻根算法的问题。我想我可以在一个单独的问题中问这个问题。

谢谢...

Answer 1

我建议使用求解器是最好的主意，主要是因为考虑各种方程的各种稳定性和收敛区域可能很棘手，而且重新发明轮子是没有用的。虽然我从未真正解决过像您提到的那样的系统，但我认为以下库之一很可能有您想要的解决方案:

• http://www.gnu.org/software/gsl/
• http://www.scimath.com/
• http://www.moshier.net/#Cephes
• 当然还有更多内容可以通过 google 找到:) [这些是我所熟悉的]

此外，如果这些都没有完全满足您的要求，您可以了解 GNU Octave 或类似的东西如何解决该系统，然后阅读他们关于用于实现解决该系统所需函数的算法的文档。从这里开始，更多的是弄清楚如何使用算法，如何实现它，以及在哪些情况下值得这样做(Octave、Matlab、Mathematica 的文档非常全面，并列出了在大多数情况下定义了该算法的出版物) ，如果您正在寻找开源/免费替代品，还有 Scilab 和 SageMath，并且有多种方法可以在 C++ 中使用这些例程(但我不确定这有多容易或多困难)

希望有帮助。