log_add_impl.hpp

Go to the documentation of this file.

#ifndef MLPACK_CORE_MATH_LOG_ADD_IMPL_HPP
#define MLPACK_CORE_MATH_LOG_ADD_IMPL_HPP

#include "log_add.hpp"

namespace mlpack {
namespace math {

template<typename T>
T LogAdd(T x, T y)
{
  T d, r;
  if (x > y)
  {
    d = y - x;
    r = x;
  }
  else
  {
    d = x - y;
    r = y;
  }

  if (std::isinf(d) || std::isinf(r))
    return r;

  return r + log(1 + exp(d));
}

template<typename T>
typename T::elem_type AccuLog(const T& x)
{
  typename T::elem_type maxVal = max(x);
  if (maxVal == -std::numeric_limits<typename T::elem_type>::infinity())
    return maxVal;

  return maxVal + log(sum(exp(x - maxVal)));;
}

template<typename T, bool InPlace>
void LogSumExp(const T& x, arma::Col<typename T::elem_type>& y)
{
  arma::Col<typename T::elem_type> maxs;

  if (InPlace)
  {
    // Compute the maximum in each column (treating y as a column too).
    maxs = max(max(x, 1), y);

    y = maxs + log(sum(exp(x - repmat(maxs, 1, x.n_cols)), 1) +
        exp(y - maxs));
  }
  else
  {
    // Compute the maximum element in each column.
    maxs = max(x, 1);

    y = maxs + log(sum(exp(x - repmat(maxs, 1, x.n_cols)), 1));
  }

  if (maxs.has_inf())
  {
    y.replace(-std::numeric_limits<typename T::elem_type>::quiet_NaN(),
              -std::numeric_limits<typename T::elem_type>::infinity());
  }
}

template<typename T, bool InPlace>
void LogSumExpT(const T& x, arma::Col<typename T::elem_type>& y)
{
  arma::Row<typename T::elem_type> maxs;

  if (InPlace)
  {
    // Compute the maximum element in each column.
    maxs = max(max(x, 0), y.t());

    y = maxs.t() + log(sum(exp(x - repmat(maxs, x.n_rows, 1)), 0) +
        exp(y.t() - maxs)).t();
  }
  else
  {
    // Compute the maximum element in each column.
    arma::Row<typename T::elem_type> maxs = max(x, 0);

    y = (maxs + log(sum(exp(x - repmat(maxs, x.n_rows, 1)), 0))).t();
  }

  if (maxs.has_inf())
  {
    y.replace(-std::numeric_limits<typename T::elem_type>::quiet_NaN(),
              -std::numeric_limits<typename T::elem_type>::infinity());
  }
}

} // namespace math
} // namespace mlpack

#endif