Loading...
Searching...
No Matches
logspace.hpp
Go to the documentation of this file.
// --------------------------------------------------------------------------------------------------
// Copyright (c) 2006-2023, Knut Reinert & Freie Universität Berlin
// Copyright (c) 2016-2023, Knut Reinert & MPI für molekulare Genetik
// This file may be used, modified and/or redistributed under the terms of the 3-clause BSD-License
// shipped with this file and also available at: https://github.com/seqan/raptor/blob/main/LICENSE.md
// --------------------------------------------------------------------------------------------------
#pragma once
#include <cmath>
namespace raptor::logspace
{
constexpr double ln_2{0.693147180559945309417232121458176568L};
constexpr double negative_inf{-std::numeric_limits<double>::infinity()};
// Produces correct result if either term is -inf.
// Needs a check for the case where both terms are -inf.
[[nodiscard]] inline double add(double const log_x, double const log_y) noexcept
{
double const max{std::max(log_x, log_y)};
return max == negative_inf ? negative_inf : max + std::log1p(std::exp(-std::abs(log_x - log_y)));
}
template <typename... types>
[[nodiscard]] double add(double const log_x, double const log_y, types... logs) noexcept
{
return add(add(log_y, log_x), logs...);
}
// expm1 is more accurate than using exp if the difference is close to 0.
[[nodiscard]] inline double substract(double const log_x, double const log_y) noexcept
{
double const difference{log_y - log_x};
return log_x + difference > -ln_2 ? std::log(std::expm1(difference)) : std::log1p(-std::exp(difference));
}
struct add_fn
{
[[nodiscard]] double operator()(double const log_x, double const log_y) const noexcept
{
return add(log_x, log_y);
};
};
} // namespace raptor::logspace
cmath
std::exp
T exp(T... args)
std::expm1
T expm1(T... args)
std::numeric_limits::infinity
T infinity(T... args)
std::log1p
T log1p(T... args)
std::log
T log(T... args)
raptor::logspace::substract
double substract(double const log_x, double const log_y) noexcept
The log of a difference of two log terms. (log_x - log_y).
Definition logspace.hpp:41
raptor::logspace::add
double add(double const log_x, double const log_y) noexcept
The log of a sum of two log terms.
Definition logspace.hpp:26
std::max
T max(T... args)
raptor::logspace::add_fn
Definition logspace.hpp:48