// graph-tool -- a general graph modification and manipulation thingy // // Copyright (C) 2006-2020 Tiago de Paula Peixoto // // This program is free software; you can redistribute it and/or modify it under // the terms of the GNU Lesser General Public License as published by the Free // Software Foundation; either version 3 of the License, or (at your option) any // later version. // // This program is distributed in the hope that it will be useful, but WITHOUT // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more // details. // // You should have received a copy of the GNU Lesser General Public License // along with this program. If not, see . #ifndef UTIL_HH #define UTIL_HH #include "config.h" #include #include "cache.hh" namespace graph_tool { using namespace boost; template [[gnu::const]] inline double lbinom(T1 N, T2 k) { if (N == 0 || k == 0 || k >= N) return 0; assert(N > 0); assert(k > 0); return ((std::lgamma(N + 1) - std::lgamma(k + 1)) - std::lgamma(N - k + 1)); } template [[gnu::const]] inline double lbinom_fast(T1 N, T2 k) { if (N == 0 || k == 0 || k > N) return 0; return ((lgamma_fast(N + 1) - lgamma_fast(k + 1)) - lgamma_fast(N - k + 1)); } template [[gnu::const]] inline double lbinom_careful(T1 N, T2 k) { if (N == 0 || k == 0 || k >= N) return 0; double lgN = std::lgamma(N + 1); double lgk = std::lgamma(k + 1); if (lgN - lgk > 1e8) { // We have N >> k. Use Stirling's approximation: ln N! ~ N ln N - N // and reorder return - N * log1p(-k / N) - k * log1p(-k / N) - k - lgk + k * log(N); } else { return lgN - std::lgamma(N - k + 1) - lgk; } } template [[gnu::const]] inline auto lbeta(T x, T y) { return (std::lgamma(x) + std::lgamma(y)) - std::lgamma(x + y); } template T log_sum(T a, T b) { if (a < b) std::swap(a, b); return a + std::log1p(exp(b-a)); } template void remove_element(Vec& vec, PosMap& pos, Val&& val) { auto& back = vec.back(); auto& j = pos[back]; auto i = pos[val]; vec[i] = back; j = i; vec.pop_back(); } template void add_element(Vec& vec, PosMap& pos, Val&& val) { pos[val] = vec.size(); vec.push_back(val); } template bool has_element(Vec& vec, PosMap& pos, Val&& val) { size_t i = pos[val]; return (i < vec.size() && vec[i] == val); } } // namespace graph_tool #endif // UTIL_HH