階乗テーブル
(library/math/factorial.hpp)

• View this file on GitHub
• Last update: 2024-01-30 20:57:42+09:00
• Include: #include "library/math/factorial.hpp"

階乗テーブル

Required by

• Binomial Coefficient Sum (library/math/binomial_coefficient_sum.hpp)
• Factorial Large (library/math/factorial_large.hpp)
• $\displaystyle \sum _ i i ^ d r ^ i$ (library/math/sum_i^d_r^i.hpp)
• Util (library/math/util.hpp)
• Compose Exp (library/polynomial/compose_exp.hpp)
• Polynomial Taylor Shift (library/polynomial/polynomial_taylor_shift.hpp)
• library/polynomial/rook_polynomial.hpp
• Shift of Sampling Points of Polynomial (library/polynomial/shift_of_sampling_points.hpp)
• Bell Number (library/sequence/bell_number.hpp)
• Bernoulli Number (library/sequence/bernoulli_number.hpp)
• Binomial Coefficient Small Prime Mod (library/sequence/binomial_coefficient_small_prime_mod.hpp)
• Eulerian Number (library/sequence/eulerian_number.hpp)
• Stirling Number1 (library/sequence/stirling_number1.hpp)
• Stirling Number of the First Kind (Small Prime Mod) (library/sequence/stirling_number1_small_prime_mod.hpp)
• Stirling Number2 (library/sequence/stirling_number2.hpp)
• Stirling Number2 Small Prime Mod (library/sequence/stirling_number2_small_prime_mod.hpp)
• library/sequence/sum_of_powers.hpp
• test/src/math/factorial_large/yuki502.cpp

Verified with

• test/src/linear_algebra/circulant_matrix/arc139_e.test.cpp
• test/src/math/array_set_power_series/abc253_h.test.cpp
• test/src/math/array_set_power_series/abc253_h_2.test.cpp
• test/src/math/binomial_coefficient_sum/yuki2512.test.cpp
• test/src/math/factorial/binomial_coefficient_prime_mod.test.cpp
• test/src/math/factorial_large/dummy.test.cpp
• test/src/math/factorial_large/factorial.test.cpp
• test/src/math/factorial_large/many_factorials.test.cpp
• test/src/math/set_power_series/abc253_h.test.cpp
• test/src/math/set_power_series/abc253_h_2.test.cpp
• test/src/math/sum_i^d_r^i/sum_of_exponential_times_polynomial.test.cpp
• test/src/math/sum_i^d_r^i/sum_of_exponential_times_polynomial_limit.test.cpp
• test/src/math/util/abc240_g.test.cpp
• test/src/polynomial/compose_exp/arc154_f.test.cpp
• test/src/polynomial/polynomial_taylor_shift/polynomial_taylor_shift.test.cpp
• test/src/polynomial/polynomial_taylor_shift/polynomial_taylor_shift_2.test.cpp
• test/src/polynomial/prod_f(r^k_x)/yuki2097.test.cpp
• test/src/polynomial/rook_polynomial/abc272_h.test.cpp
• test/src/polynomial/rook_polynomial/dummy.test.cpp
• test/src/polynomial/shift_of_sampling_points/shift_of_sampling_points_of_polynomial.test.cpp
• test/src/sequence/bernoulli_number/bernoulli_number.test.cpp
• test/src/sequence/bernoulli_number/bernoulli_number_2.test.cpp
• test/src/sequence/eulerian_number/yuki2005-2-2.test.cpp
• test/src/sequence/eulerian_number/yuki2005-2.test.cpp
• test/src/sequence/eulerian_number/yuki2005.test.cpp
• test/src/sequence/stirling_number1/abc247_h.test.cpp
• test/src/sequence/stirling_number1/stirling_number1.test.cpp
• test/src/sequence/stirling_number1/stirling_number1_2.test.cpp
• test/src/sequence/stirling_number1_small_prime_mod/stirling_number_of_the_first_kind_small_p_large_n.test.cpp
• test/src/sequence/stirling_number2/stirling_number2.test.cpp
• test/src/sequence/stirling_number2/stirling_number2_2.test.cpp
• test/src/sequence/stirling_number2_small_prime_mod/stirling_number_of_the_second_kind_small_p_large_n.test.cpp

Code

#ifndef SUISEN_FACTORIAL
#define SUISEN_FACTORIAL

#include <cassert>
#include <vector>

namespace suisen {
    template <typename T, typename U = T>
    struct factorial {
        factorial() = default;
        factorial(int n) { ensure(n); }

        static void ensure(const int n) {
            int sz = _fac.size();
            if (n + 1 <= sz) return;
            int new_size = std::max(n + 1, sz * 2);
            _fac.resize(new_size), _fac_inv.resize(new_size);
            for (int i = sz; i < new_size; ++i) _fac[i] = _fac[i - 1] * i;
            _fac_inv[new_size - 1] = U(1) / _fac[new_size - 1];
            for (int i = new_size - 1; i > sz; --i) _fac_inv[i - 1] = _fac_inv[i] * i;
        }

        T fac(const int i) {
            ensure(i);
            return _fac[i];
        }
        T operator()(int i) {
            return fac(i);
        }
        U fac_inv(const int i) {
            ensure(i);
            return _fac_inv[i];
        }
        U binom(const int n, const int r) {
            if (n < 0 or r < 0 or n < r) return 0;
            ensure(n);
            return _fac[n] * _fac_inv[r] * _fac_inv[n - r];
        }
        template <typename ...Ds, std::enable_if_t<std::conjunction_v<std::is_integral<Ds>...>, std::nullptr_t> = nullptr>
        U polynom(const int n, const Ds& ...ds) {
            if (n < 0) return 0;
            ensure(n);
            int sumd = 0;
            U res = _fac[n];
            for (int d : { ds... }) {
                if (d < 0 or d > n) return 0;
                sumd += d;
                res *= _fac_inv[d];
            }
            if (sumd > n) return 0;
            res *= _fac_inv[n - sumd];
            return res;
        }
        U perm(const int n, const int r) {
            if (n < 0 or r < 0 or n < r) return 0;
            ensure(n);
            return _fac[n] * _fac_inv[n - r];
        }
    private:
        static std::vector<T> _fac;
        static std::vector<U> _fac_inv;
    };
    template <typename T, typename U>
    std::vector<T> factorial<T, U>::_fac{ 1 };
    template <typename T, typename U>
    std::vector<U> factorial<T, U>::_fac_inv{ 1 };
} // namespace suisen

#endif // SUISEN_FACTORIAL

#line 1 "library/math/factorial.hpp"



#include <cassert>
#include <vector>

namespace suisen {
    template <typename T, typename U = T>
    struct factorial {
        factorial() = default;
        factorial(int n) { ensure(n); }

        static void ensure(const int n) {
            int sz = _fac.size();
            if (n + 1 <= sz) return;
            int new_size = std::max(n + 1, sz * 2);
            _fac.resize(new_size), _fac_inv.resize(new_size);
            for (int i = sz; i < new_size; ++i) _fac[i] = _fac[i - 1] * i;
            _fac_inv[new_size - 1] = U(1) / _fac[new_size - 1];
            for (int i = new_size - 1; i > sz; --i) _fac_inv[i - 1] = _fac_inv[i] * i;
        }

        T fac(const int i) {
            ensure(i);
            return _fac[i];
        }
        T operator()(int i) {
            return fac(i);
        }
        U fac_inv(const int i) {
            ensure(i);
            return _fac_inv[i];
        }
        U binom(const int n, const int r) {
            if (n < 0 or r < 0 or n < r) return 0;
            ensure(n);
            return _fac[n] * _fac_inv[r] * _fac_inv[n - r];
        }
        template <typename ...Ds, std::enable_if_t<std::conjunction_v<std::is_integral<Ds>...>, std::nullptr_t> = nullptr>
        U polynom(const int n, const Ds& ...ds) {
            if (n < 0) return 0;
            ensure(n);
            int sumd = 0;
            U res = _fac[n];
            for (int d : { ds... }) {
                if (d < 0 or d > n) return 0;
                sumd += d;
                res *= _fac_inv[d];
            }
            if (sumd > n) return 0;
            res *= _fac_inv[n - sumd];
            return res;
        }
        U perm(const int n, const int r) {
            if (n < 0 or r < 0 or n < r) return 0;
            ensure(n);
            return _fac[n] * _fac_inv[n - r];
        }
    private:
        static std::vector<T> _fac;
        static std::vector<U> _fac_inv;
    };
    template <typename T, typename U>
    std::vector<T> factorial<T, U>::_fac{ 1 };
    template <typename T, typename U>
    std::vector<U> factorial<T, U>::_fac_inv{ 1 };
} // namespace suisen

Back to top page