二項係数の prefix sum の内部演算（高速版）

components/binomial_prefix_core.hpp · GitHub で見る

全体サイト（256 KiB向け） · 高速実装

軽量実装（64 KiB）はこちら

二項係数の prefix sum の共通部品です。binomial_prefix_detail::Table<Mint> を使います。

操作・メンバ	内容・返り値	計算量
Table(n=0)	\(0\leq i\leq n\) の階乗・逆階乗・2 の冪を前計算	\(O(n)\) と逆元 1 回
build(n)	前計算を作り直す	\(O(n)\) と逆元 1 回
binom(n,k)	\(\displaystyle \binom nk\)。k < 0 または k > n なら 0	\(O(1)\)
center(n)	\(\displaystyle \sum_{i=0}^{\lfloor(n+1)/2\rfloor-1}\binom ni\)	\(O(1)\)
range(n,l,r)	\(\displaystyle \sum_{i=l}^{r-1}\binom ni\)。0 <= l <= r <= n+1 が必要	\(O(r-l)\)
factorial	階乗の vector<Mint>	\(O(1)\)
inverse_factorial	逆階乗の vector<Mint>	\(O(1)\)
power2	2 の冪の vector<Mint>	\(O(1)\)
half	奇素数 mod では 2 の逆元	\(O(1)\)

法は素数、前計算の上限は非負で INT_MAX 未満かつ法未満にします。各操作の n は前計算の範囲内で指定します。Mint の条件は公開ライブラリと同じです。

Mint = mint では modint の階乗・逆階乗生成を使い、他の型ではスカラー処理を使います。NTT・FPS には依存しません。

前計算は \(O(n)\) 時間・空間と逆元計算 1 回です。binom・center は \(O(1)\)、range は \(O(r-l)\) 時間・\(O(1)\) 追加空間です。

依存ライブラリ

• modint998244353 の共通剰余演算
• modint998244353 の階乗・逆階乗生成
• modint998244353

検証コード

• test/binomial_prefix_sum.test.cpp
• test/binomial_prefix_sum_online.test.cpp

コード

#ifndef BINOMIAL_PREFIX_CORE_HPP
#define BINOMIAL_PREFIX_CORE_HPP

#include "../math/mint.hpp"
#include "mint_factorials.hpp"
#include <cassert>
#include <cstdint>
#include <limits>
#include <type_traits>
#include <vector>

namespace binomial_prefix_detail {
template<class Mint> struct Table {
	std::vector<Mint> factorial, inverse_factorial, power2;
	Mint half;
	explicit Table(int n = 0) { build(n); }
	void build(int n) {
		assert(0 <= n && n < std::numeric_limits<int>::max() && std::uint64_t(n) < Mint::mod() && Mint::mod() > 1);
		factorial.resize(n + 1), inverse_factorial.resize(n + 1), power2.resize(n + 1);
		if constexpr (std::is_same_v<Mint, mint>) {
			mint_detail::factorials(factorial.data(), inverse_factorial.data(), n + 1);
			power2[0] = Mint(1);
			for (int i = 1; i <= n; ++i) power2[i] = power2[i - 1] + power2[i - 1];
		} else {
			factorial[0] = power2[0] = Mint(1);
			for (int i = 1; i <= n; ++i) factorial[i] = factorial[i - 1] * Mint(i), power2[i] = power2[i - 1] + power2[i - 1];
			inverse_factorial[n] = factorial[n].inv();
			for (int i = n; i > 0; --i) inverse_factorial[i - 1] = inverse_factorial[i] * Mint(i);
		}
		half = Mint((std::uint64_t(Mint::mod()) + 1) / 2);
	}
	Mint binom(int n, int k) const {
		if (k < 0 || k > n) return Mint(0);
		return factorial[n] * inverse_factorial[k] * inverse_factorial[n - k];
	}
	Mint center(int n) const {
		if (!n) return Mint(0);
		return n & 1 ? power2[n - 1] : (power2[n] - binom(n, n / 2)) * half;
	}
	Mint range(int n, int l, int r) const {
		assert(0 <= l && l <= r && r <= n + 1);
		if constexpr (std::is_same_v<Mint, mint>) {
			__uint128_t sum = 0;
			for (int i = l; i < r; ++i) sum += std::uint64_t(inverse_factorial[i].v) * inverse_factorial[n - i].v;
			return factorial[n] * Mint::word(Mint::reduce(std::uint64_t(sum % Mint::mod())));
		} else {
			Mint sum(0);
			for (int i = l; i < r; ++i) sum += inverse_factorial[i] * inverse_factorial[n - i];
			return factorial[n] * sum;
		}
	}
};
}

#endif

Back to top page