#ifndef log_add_hh #define log_add_hh #include <limits> #include <iostream> #include <cassert> #include <cmath> template <typename T> struct Log { static T zero() { return -std::numeric_limits<T>::infinity(); } static T add(T l1, T l2) { if (l1 == zero()) return l2; if (l1 > l2) return l1 + std::log(1 + exp(l2 - l1)); else return l2 + std::log(1 + exp(l1 - l2)); } static T subtract(T l1, T l2) { //std::assert(l1 >= l2); return l1 + log(1 - exp(l2 - l1)); } }; #endif