c - C 中的 logsumexp 实现？

Question

有人知道提供 logsumexp 函数的开源数值 C 库吗？

logsumexp(a) 函数计算数组 a 的分量的指数 log(e^{a_1}+...e^{a_n}) 之和，避免数值溢出。

Answer 1

这是一个非常简单的从头开始的实现(至少经过最低限度的测试):

double logsumexp(double nums[], size_t ct) {
  double max_exp = nums[0], sum = 0.0;
  size_t i;

  for (i = 1 ; i < ct ; i++)
    if (nums[i] > max_exp)
      max_exp = nums[i];

  for (i = 0; i < ct ; i++)
    sum += exp(nums[i] - max_exp);

  return log(sum) + max_exp;
}

这是有效地将所有参数除以最大参数的技巧，然后在末尾添加其对数以避免溢出，因此添加大量类似缩放的值时表现良好，但错误会逐渐增加如果某些参数比其他参数大很多个数量级。

如果你想让它在给定 0 个参数时运行而不会崩溃，你必须为此添加一个案例:)