nCr mod p ($r \leq 10^7$, クエリ $O(1)$ or $O(r)$) (combinatorics/ncr2.cpp)

• View this file on GitHub
• View document part on GitHub
• Last update: 2023-10-06 03:54:40+09:00

なにこれ

${}_nC_r\ (\mathrm{mod}\ p)$ を求める． 計算量が $r$ 依存．

コンストラクタ

• NCR(n, r, p)：$r$ の最大値を r，法を p とする．($n$ は固定値)
• NCR(r, p)：$r$ の最大値を r，法を p とする．

$n \leq 10^9$，$r \leq 10^7$ 程度．$p$ は素数．

メンバ関数

• calc(n, r)：計算量 $O(1)$ or $O(r)$
  ${}_nC_r$ を求める．
  $n$ が固定値のときは $O(1)$，不定値のときは $O(r)$

計算量

• NCR(n, r, p)：$O(r)$
• NCR(r, p)：$O(r)$
• calc(n, r)：$n$ の最大値が既知の場合 $O(1)$，未知の場合 $O(r)$

Verified with

• test/aoj_DPL_5_E.cpp
• test/aoj_DPL_5_E_noset_n.cpp

Code

#pragma once

#ifndef call_include
#define call_include
#include <bits/stdc++.h>
using namespace std;
#endif

struct NCR {
private:
	vector<long long> comb, inv, finv;
	long long P;

	void calc_inv(long long r) {
		inv[1] = finv[0] = finv[1] = 1LL;
		for(int i = 2; i <= r; i++) {
			inv[i] = P - inv[P % i] * (P / i) % P;
			finv[i] = finv[i - 1] * inv[i] % P;
		}
	}

public:
	NCR(long long n, long long r, long long p) : comb(n + 1), inv(r + 1), finv(r + 1), P(p) {
		if(n / 2 < r) r = n / 2;
		calc_inv(r);

		comb[0] = 1;
		for(int i = 1; i <= r; i++) {
			comb[i] = comb[i - 1] * (n - i + 1) % P * inv[i] % P;
		}
	}

	NCR(long long r, long long p) : inv(r + 1), finv(r + 1), P(p) {
		calc_inv(r);
	}

	long long calc(long long n, long long r) {
		assert(r >= 0);
		if(r > n) return 0;
		if(r > n / 2) r = n - r;
		if(comb.size() > 0) return comb[r];
		else {
			long long f = 1;
			for(long long i = n; i > n - r; i--) (f *= i) %= P;
			return f * finv[r] % P;
		}
	}
};

#line 2 "combinatorics/ncr2.cpp"

#ifndef call_include
#define call_include
#include <bits/stdc++.h>
using namespace std;
#endif

struct NCR {
private:
	vector<long long> comb, inv, finv;
	long long P;

	void calc_inv(long long r) {
		inv[1] = finv[0] = finv[1] = 1LL;
		for(int i = 2; i <= r; i++) {
			inv[i] = P - inv[P % i] * (P / i) % P;
			finv[i] = finv[i - 1] * inv[i] % P;
		}
	}

public:
	NCR(long long n, long long r, long long p) : comb(n + 1), inv(r + 1), finv(r + 1), P(p) {
		if(n / 2 < r) r = n / 2;
		calc_inv(r);

		comb[0] = 1;
		for(int i = 1; i <= r; i++) {
			comb[i] = comb[i - 1] * (n - i + 1) % P * inv[i] % P;
		}
	}

	NCR(long long r, long long p) : inv(r + 1), finv(r + 1), P(p) {
		calc_inv(r);
	}

	long long calc(long long n, long long r) {
		assert(r >= 0);
		if(r > n) return 0;
		if(r > n / 2) r = n - r;
		if(comb.size() > 0) return comb[r];
		else {
			long long f = 1;
			for(long long i = n; i > n - r; i--) (f *= i) %= P;
			return f * finv[r] % P;
		}
	}
};

Back to top page